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Abstract 

04 \ Given a point s and a set of h pairwise disjoint polygonal obstacles of totally n vertices in 

■ the plane, we present a new algorithm for building an Li shortest path map of size 0(n) in 0{T) 
\ time and 0{n) space such that for any query point t, the length of the Li shortest obstacle- 
avoiding path from s to t can be reported in 0(log n) time and the actual shortest path can be 

I found in additional time proportional to the number of edges of the path, where T is the time for 

. triangulating the free space. It is currently known that T — 0{n + h\o^^^ h) for an arbitrarily 

small constant e > 0. If the triangulation can be done optimally (i.e., T = 0{n + ft, log /i)), then 
our algorithm is optimal. Previously, the best algorithm computes such an Li shortest path map 
in 0(rt log n) time and 0{n) space. Our techniques can be extended to obtain improved results 
for other related problems, e.g., computing the Li geodesic Voronoi diagram for a set of point 

■ sites in a polygonal domain, finding shortest paths with fixed orientations, finding approximate 
. Euclidean shortest paths, etc. 

O . 

1 Introduction 

^ ■ Computing obstacle-avoiding shortest paths in the plane is a fundamental problem in computational 

■ geometry and has many applications. The Euclidean version that measures the path length by the 
Euclidean distance has been well studied (e.g., see [7, 8, 15, 19, 22, 24, 25, 31, 33, 35]). In this 

\ paper, we consider the Li version, defined as follows. Given a point s and a set of h pairwise disjoint 

CnJ ■ polygonal obstacles, V = {Pi, • • • > Ph}-, of totally n vertices in the plane, where s is considered 

■ as a special point obstacle, the plane minus the interior of the obstacles is called the free space of 
V. Two obstacles are pairwise disjoint if they do not intersect in their interior. The Li shortest 
path map problem, denoted by Li-SPM, is to compute a single-source shortest path map (SPM for 
short) with s as the source point such that for any query point t, an Li shortest obstacle-avoiding 
path from s to t can be obtained efficiently. Note that such a path can consist of any polygonal 

. segments but the length of each segment of the path is measured by the Li metric. 

We say that an SPM has standard query performances if for any query point t, the length of 
the Li shortest obstacle-avoiding path from s to t can be reported in O(logn) time and an actual 
shortest path can be found in additional time proportional to the number of edges (or turns) of the 
path. 

If the input also includes another point t and the problem only asks for one single Li shortest 
path from s to t, then we call this problem version the Li shortest path problem, denoted by Li-SP. 

A closely related problem version solvable by our approach is to find shortest rectilinear paths. 
A rectilinear path is a path each of whose edges is parallel to a coordinate axis and its length 
is measured by the Euclidean distances or Li distances of its segments (they are the same for 
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rectilinear paths). Rectilinear shortest paths are used widely in VLSI design and network wire- 
routing applications. As shown in [10, 27, 29, 30], it is easy to convert an arbitrary polygonal 
path to a rectilinear path with the same Li length. Thus, in this paper, we focus on computing 
polygonal paths measured by the Li distance. 

1.1 Previous Work 

The Li-SP problem has been studied extensively (e.g., see [6, 10, 11, 27, 29, 30, 36]). In general, 
there are two approaches for solving this problem: Constructing a sparse "path preserving" graph 
(analogous to a visibility graph), and the continuous Dijkstra paradigm. Clarkson, Kapoor, and 
Vaidya [10] constructed a graph of 0(n log n) nodes and 0(n log n) edges such that a shortest path 
can be found in the graph in 0{n log^ n) time; subsequently, they gave an algorithm of 0{n log^"'' n) 
time and 0(n log^'^ n) space [11]. Based on some observations, Chen, Klenk, and Tu [6] showed that 
the problem was solvable in 0(n log^'^ ra) time and 0(n log n) space. By applying the continuous 
Dijkstra paradigm, Mitchell [29, 30] solved the problem in O(nlogn) time and 0{n) space. An 
0{n + hlogh) time lower bound can be established for solving Li-SP (e.g., based on the results 
in [12]). Hence, Mitchell's algorithm is worst-case optimal. Recently, by using a corridor structure 
and building a smaller path preserving graph, Inkulu and Kapoor [21] solved the Li-SP problem 
in 0{n + /ilog^'^ n) time and 0{n + hlog^'^ h) space. 

For the query version of the problem, i.e., Li-SPM, Mitchell's algorithm [29, 30] builds an SPM 
of size 0{n) in 0(n log n) time and 0{n) space with the standard query performances. 

In addition, for the convex case where all polygonal obstacles in V are convex, to our best 
knowledge, we are not aware of any previous better results than those mentioned above. 

1.2 Our Results 

We present an algorithm for Li-SPM that builds an SPM of size 0{n) in 0{T) time and 0{n) 
space with the standard query performances, where T always refers to the time for triangulating 

the free space of V in the paper. It is obvious to sec that given an SPM, we can always add h ~ 1 
line segments in the free space to connect the obstacles in V together to obtain a single simple 
polygon and then triangulate the free space, in totally 0{n) time [2, 3]. It is currently known 
that T = Q(n + hlogh) and T = 0{n + hlog^'^^ h) [2], where e is an arbitrarily small positive 
constant. Therefore, we essentially solve Li-SPM in 0(T) time. In other words, our result shows 
that building an SPM is equivalent to triangulating the free space of V in terms of the running 
time. 

Our approach uses Mitchell's algorithm [29, 30] as a procedure and further explores the corridor 
structure of P [25]. One interesting observation we found is that to find an Li shortest path among 
convex obstacles, it is sufficient to consider only the at most four extreme vertices (along the 
horizontal and vertical directions) of each obstacle (these vertices define a core for each obstacle). 
Mitchell's algorithm is then applied to these cores, which takes only 0{hlogh) time. More work 
needs to be done for computing an SPM. For example, one key result we have is that we give an 
0{n' + m!) time algorithm for a special case of constructing the Li geodesic Voronoi diagram in a 
simple polygon of n' vertices for m! weighted point sites, where the sites all lie outside the polygon 
and influence the polygon through one (open) edge (see Fig. 1). We are not aware of any specific 
previous work on this problem, although an 0{{n' +m') \og{n' +m')) time solution may be obtained 
by standard techniques. Our linear time algorithm, which is clearly optimal, may be interesting in 
its own right. 
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Figure 1: (a) Three weighted sites (in red) and a simple polygon B with an open edge cd. The goal is to compute 
the Li geodesic Voronoi diagram in B with respect to the three sites which influence B only through the edge cd. 
(b) Illustrating a possible solution: B is partitioned into three Voronoi regions VD{ri) for each r^, 1 < i < 3. 

For the convex case where all obstacles in V are convex, we can find a shortest s-t path in 
0{n + hlogh) time and 0{n) space since the triangulation can be done in 0{n + hlogh) time (e.g., 
by the approaches in [2, 20]); this is optimal. A by-product of our techniques, which may be a 
little "surprising", is that in 0{n-\- hlogh) time and 0{n) space, we can build an SPM of size 0{h) 
(instead of 0{n)) such that the shortest path length queries are answered in 0{logh) time each 
(instead of O(logn) time). 

1.3 Applications 

Our techniques can be extended to solve other problems. 

The Li geodesic Voronoi diagram problem, denoted by Li-GVD, is defined as follows. Given an 
obstacle set V and a set of m point sites in the free space, compute the geodesic Voronoi diagram 
for the m point sites under the Li distance metric among the obstacles in V. Mitchell [29, 30], 
solves the Li-GVD problem in 0((n + m) log(n + m)) time. Our approach can compute it in 
0{T' + n + (m + h) log(m + h)) time, where T' is the time for triangulating the free space along 
with the m point sites. It is known that T' = 0{n -\- (m + h) log"^^'^(m + h)) [2] or alternatively we 
can obtain T' = 0{n + hlog^'^'' h + mlogn). Note that when applying our algorithm to a single 
simple polygon P of n vertices, the Li geodesic Voronoi diagram for m point sites in P can be 
obtained in 0{n + mlog^''''^m) or 0{n + m(logn + logm)) time. In comparison, the Euclidean 
version of the one simple polygon case was solved in 0{{n + m) log(n + m)) time [32]. 

We also give better results for the shortest path problem in "fixed orientation metrics" [29, 30, 
37], for which a sought path is allowed to follow only a given set of orientations. For a number c 
of given orientations, Mitchell's algorithm [29, 30] finds such a shortest path in 0(cn log n) time 
and 0{cn) space, and our algorithm takes 0{n + hlog^^"^ h + c^hlogch) time and 0{n + (?h) 
space. In addition, our approach also leads to an 0{n + hlog}^'^ h + {l/5)h\og time algorithm 
for computing a (5-optimal Euclidean shortest path among polygonal obstacles for any constant 
6 > 0. For this problem, Mitchell's algorithm [29, 30] takes 0{{^/T/5)n log n) time, and Clarkson's 
algorithm [9] runs in 0{{l/6)n\ogn) time. 

2 An Overview of Our Approaches 

In this section, we give an overview of our approaches as well as the organization of this paper. 
Denote by the free space of V. We begin with our algorithm for the convex case, which is a key 
procedure for solving the general problem. 
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We first discuss tlic Li-SP problem. In the convex case, each obstacle in V = {Pi,P2, ■ ■ ■ ,Ph} 
is convex. For each Pi £ V, we compute its core, denoted by core{Pi), which is a simple polygon 
by connecting the topmost, leftmost, bottommost, and rightmost points of Pi. Let core{V) be 
the set of all h cores of V. For any point t in the free space T, we show that given any shortest 
s-t path avoiding all cores in core{V), wc can find in 0(n) time a shortest s-t path avoiding all 
obstacles in V with the same Li length. Based on this observation, our algorithm has two main 
steps: (1) Apply Mitchell's algorithm [29, 30] on core{V) to compute a shortest s-t path -Kcoreis^t) 
avoiding the cores in core{V), which takes 0{hlogh) time since each core in core{V) has at most 
four vertices; (2) based on TTcoreis,t), compute a shortest s-t path avoiding all obstacles in V in 
0(n) time. This algorithm takes overall 0(n + hlogh) time and 0(n) space. 

To build an SPM in T (with respect to the source point s), similarly, we first apply Mitchell's 
algorithm on core(V) to compute an SPM of 0{h) size in the free space with respect to all cores, 
which can be done in 0(n -\- hlogh) time and 0(n) space. Based on the above SPM, in additional 
0{n) time, we are able to compute an SPM in F. Our results for the convex case are given in 
Section 3. 

For the general problem where the obstacles in V are not necessarily convex, based on a tri- 
angulation of the free space JF, we first compute a corridor structure [25], which consists of 0{h) 
corridors and 0{h) junction triangles. Each corridor possibly has a corridor path. As in [25], the 
corridor structure can be used to partition the plane into a set V' of 0{h) pairwise disjoint convex 
polygons of totally 0(n) vertices such that a shortest s-t path in is a shortest s-t path avoiding 
the convex polygons in V' and possibly containing some corridor paths. All corridor paths are 
contained in the polygons of V' . Thus, in addition to the corridor paths, finding a shortest path 
is reduced to an instance of the convex case. By incorporating the corridor path information into 
Mitchell's continuous Dijkstra paradigm [29, 30], our algorithm for the convex case can be modified 
to find a shortest path in 0{T) time. The above algorithm is presented in Section 4. 

Sections 4.3, 5, and 6 are together devoted to compute an SPM in F (Section 4.3 outlines the 
algorithm) . We use the corridor structure to partition T into the ocean M. , bays, and canals. While 
the ocean Ai may be multiply connected, every bay or canal is a simple polygon. Each bay has a 
single common boundary edge with Ai and each canal has two common boundary edges with Ai. 
But two bays or two canals, or a bay and a canal do not share any boundary edge. A common 
boundary edge of a bay (or canal) with M is called a gate. Thus each bay has one gate and each 
canal has two gates. Further, the ocean A4 is exactly the free space with respect to the convex 
polygonal set V' . By modifying our algorithm for the convex case, we can compute an SPM in A4 
in 0{T) time. This part is discussed in Section 4.3. 

Denote by SPM{M) the SPM in M. To obtain an SPM in J", we need to "expand" SPM{M) 
into all bays and canals through their gates. Here, a challenging subproblem is to solve efficiently a 
special case of the (additively) weighted Li geodesic Voronoi diagram problem on a simple polygon 
B: The weighted point sites all lie outside B and influence B through one (open) edge (e.g., see 
Fig. 1). The subproblem models the procedure of expanding SPM{M.) into a bay, where the polygon 
B is the bay, the point sites are obstacle vertices in Ai, the weight of each site is the length of its 
shortest path to the source point s, and the edge of the polygon (e.g., cd in Fig. 1) is the gate of 
the bay. As discussed before, we give a linear time solution for this subproblem in Section 5. Note 
that although our presentation for solving the subproblem is long and technically complicated, the 
algorithm itself is simple and easy to implement; our effort is mostly for simplifying the algorithm 
and showing its correctness. 

Expanding SPM{Ai) into canals, which is discussed in Section 6, is also done in linear time by 
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using our solution for the above subproblem as a main procedure. In summary, given SPM(Ai), 
computing an SPM for the entire free space takes additional 0(n) time. 

We discuss a little more about the above challenging subproblem. The problem may not look 
"challenging" at all as it can be solved by many existing techniques. For example, one may attempt 
to use the continuous Dijkstra approach [29, 30] to let the "wavelet" enter into the bays/canals. 
However, that would lead to an 0{{n' + m') log(n' + m')) time solution for the subproblem since it 
takes logarithmic time to process each event, where n' is the number of vertices of B and m' is the 
number of weighted sites, and consequently it would take an overall 0(n log n) time for building 
an SPM in J^. One may also want to use a sweeping algorithm [14], which would also lead to 
an 0((n' + m') log(n' + m')) time solution since again it takes logarithmic time to process each 
event. In addition, the divide-and-conquer approach [34] would also take 0((n' + m') log(ra' + m')) 
time since the merge procedure takes linear time. Our algorithm for the subproblem, which can 
be viewed as an incremental approach, takes 0{n' + m') time. Incremental approaches have been 
widely used in geometric algorithms, and normally they can result in good randomized algorithms. 
Incremental approaches have also been used for constructing Voronoi diagrams, which usually take 
quadratic time. Our result demonstrates that incremental approaches are able to yield optimal 
deterministic solutions for building Voronoi diagrams, and the success of it hinges on discovering 
many geometric properties of the problem. We should point out that our techniques for solving the 
challenging subproblem are quite independent of other parts of the paper. 

In Section 7, we generalize our techniques to solve some related problems discussed in Section 
1.3. Section 8 concludes the paper. 

As in [29, 30], for simplicity of discussion, we assume that the free space T is connected and 
the point t is always in (thus, a feasible s-t path always exists), and no two obstacle vertices lie 
on the same horizontal or vertical line. In the rest of this paper, unless otherwise stated, a shortest 
path always refers to an Li shortest path and a length is always in the Li metric. 

3 Shortest Paths among Convex Obstacles 

In this section, we give our algorithms for the convex case, which are also used for the general 
case in later sections. Let V' = {P[,P2 ■ ■ ■ be a set of h pairwise disjoint convex polygonal 
obstacles of totally n vertices. With respect to the source point s, our algorithm builds an SPM of 
0(n) size with standard query performances in 0(n + hlogh) time and 0{n) space. 

3.1 Notation and Observations 

For each convex polygon P- G V' , we define its core, denoted by core{P-), as the simple polygon 
by connecting the leftmost, topmost, rightmost, and bottommost vertices of P^ with line segments 
(see Fig. 2). Note that core{Pl) is contained in P/ and has at most four edges. Let core{V') be the 
set of the cores of all obstacles in V' . Consider a point t in the free space T. A key observation (to 
be proved) is that a shortest s-t path avoiding the cores in coreiV') corresponds to a shortest s-t 
path avoiding the obstacles in V' with the same Li length. Note that a path avoiding the cores in 
coreiV') may intersect the interior of some obstacles in V' . 

To prove the above key observation, we first define some concepts. Consider an obstacle P/ and 
core{Pl). For each edge ab of core{Pl) with vertices a and b, if ab is not an edge of P/, then it 
divides P- into two polygons, one of them containing core{P'); we call the one that docs not contain 
core{Pl) an ear of P/ based on ab, denoted by ear{ab) (see Fig. 2). If ab is also an edge of Pj, then 
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Figure 2: Illustrating the core and ears of a convex Figure 3: The line segment cd penetrates ear{ab); cd 
obstacle; ear{ab) is indicated. intersects the obstacle path of ear{ab) at e and /. 

ear{ab) is not defined. Note that ear{ab) has only one edge bounding core{Pl), i.e., ab, which we 
cah its core edge. The other edges of ear{ab) are on the boundary of P[, which we call obstacle 
edges. There are two paths between a and b along the boundary of ear{ab): One path is the core 
edge ab and the other consists of all its obstacle edges. We call the latter path the obstacle path 
of the ear. A line segment is positive-sloped (resp., negative-sloped) if its slope is positive (resp., 
negative). An ear is positive-sloped (resp., negative-sloped) if its core edge is positive-sloped (resp., 
negative-sloped). Note that by our assumption no two obstacle vertices lie on the same horizontal 
or vertical line, and thus no ear has a horizontal or vertical core edge. A point p is higher (resp., 
lower) than another point q if the y-coordinate of p is no smaller (resp., no larger) than that of q. 
The next observation is self-evident. 

Observation 1 For any ear, its obstacle path is monotone in both the x- and y- coordinates. Specif- 
ically, consider an ear ear{ah) and suppose the vertex a is lower than the vertex b. If ear{ab) is 
positive-sloped, then the obstacle path from a to b is monotonically increasing in both the x- and 
y- coordinates; if it is negative-sloped, then the obstacle path from a to b is monotonically decreasing 
in the x-coordinates and monotonically increasing in the y- coordinates. 

For an ear ear{a,b) and a line segment cd, we say that cd penetrates ear{ab) if the following hold 
(see Fig. 3): (1) cd intersects the interior of ear{ab), (2) neither c nor d is in the interior of ear{ab), 
and (3) cd does not intersect the core edge ab at its interior. The next lemma will be useful later. 

Lemma 1 Suppose a line segment cd penetrates an ear ear{ab). If cd is positive-sloped (resp., 
negative-sloped), then ear{ab) is also positive-sloped (resp., negative-sloped). 

Proof: We only prove the case when cd is positive-sloped since the other case is similar. 

Assume to the contrary that ear{ab) is negative-sloped. Without loss of generality (WLOG), we 
assume a is lower than b. By Observation 1, the obstacle path of ear{ab) from a to 6 is monotonically 
decreasing in the x-coordinates. Thus, the rightmost point and leftmost point of ear{ab) are a and 
b, respectively. Note that ear{ab) is contained in the region between the two vertical lines passing 
through a and b. Since cd is positive-sloped and ab is negative-sloped, if cd intersects an interior 
point of ear(ab), then cd must cross ab at an interior point. But since cd penetrates ear{ab), cd 
cannot intersect any interior point of ab. Hence, we have a contradiction. The lemma thus follows. 
□ 

Clearly, if cd penetrates the ear ear{ab), then cd intersects the boundary of ear{ab) at two 
points and both points lie on the obstacle path of ear{ab) (e.g., see Fig. 3). 

Lemma 2 Suppose a line segmxmt cd penetrates an ear ear{ab). Let e and f be the two points on 
the obstacle path of ear{ab) that cd intersects. Then the Li length of the line segment ef is equal 
to that of the portion of the obstacle path of ear{ab) between e and f (see Fig. 3). 
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Proof: WLOG, suppose cd is positive-sloped and e is lower than /. By Lemma 1, ear{ab) is 
also positive-sloped. The segment ef from e to / is monotonically increasing in both the x- and 
y-coordinates. Denote by ef the portion of the obstacle path of ear{ab) between e and /. Since 
ear {ah) is positive-sloped, by Observation 1, the portion ef from e to / is monotonically increasing 
in both the x- and y-coordinates. Therefore, the Li lengths of ef and ef are equal. The lemma 
thus follows. □ 

If cd penetrates ear{ab), then by Lemma 2, we can obtain another path from c to d by replacing 
ef with the portion of the obstacle path of ear {ah) between e and / such that the new path has 
the same Li length as cd and the new path does not intersect the interior of ear {ah). 

The results in the following lemma have been proved in [29, 30]. 

Lemma 3 [29, 30] There exists a shortest s-t path in the free space such that if the path makes a 
turn at a point p, then p is an obstacle vertex. 

We call a shortest path that satisfies the property in Lemma 3 a vertex-preferred shortest path. 
Mitchell's algorithm [29, 30] can find a vertex-preferred shortest s-t path. Denote by Tri{V') a 
triangulation of the free space and the space inside all obstacles. Note that the free space can be 
triangulated in 0{n + h\og h) time [2, 20] and the space inside all obstacles can be triangulated in 
totally 0(n) time [3]. Hence, TniV') can be computed in 0{n + hlogh) time. The next lemma 
gives our key observation. 

Lemma 4 Given a vertex-preferred shortest s-t path that avoids the polygons in core{V'), we can 
find in 0{n) time a shortest s-t path with the same Li length that avoids the obstacles in V' . 

Proof: Consider a vertex-preferred shortest s-t path for core{V'), denoted by 7Tcore{s,t). Suppose 
it makes turns at pi,p2, . . . ,Pki ordered from s to t along the path, and each pi is a vertex of a core 
in core{V'). Let po = s and Pk+i = t. Then for each i = 0,1, . . . ,k, the portion of 7rcore{s, t) from 
Pi to pi+i is the line segment piPi+i, which does not intersect the interior of any core in core{V'). 
Below, we first show that we can find a path from pi to pj+i such that it avoids the obstacles in V 
and has the same Li length as PiPi+i- 

If PiPi+i does not intersect the interior of any obstacle in V' , then we are done with piPi^i. 
Otherwise, because PipJ+l avoids core{V'), it can intersects only the interior of some ears. Consider 
any such ear ear {ah). Below, we prove that PiPi+T penetrates ear{ab). 

First, we already know that pipi^i intersects the interior of ear{ab). Second, it is obvious that 
neither pi nor pi+i is in the interior of ear{ab). It remains to show that piPi+i cannot intersect the 
core edge ab of ear{ab) at the interior of ab. Denote by A' € V' the obstacle that contains ear{ab). 
The interior of ab is in the interior of A'. Since PiPi+T does not intersect the interior of A', PiPi+T 
cannot intersect ab at its interior. Therefore, PiPi+T penetrates ear{ab). 

Recall that we have assumed that no two obstacle vertices lie on the same horizontal or vertical 
line. Since both pi and pi+i are obstacle vertices, the segment PiPi+T is either positive-sloped or 
negative-sloped. WLOG, assume PipJ+l is positive-sloped. By Lemma 1, ear{ab) is also positive- 
sloped. Let e and / denote the two intersection points between PiPi+i and the obstacle path of 
ear{ab), and ef denote the portion of the obstacle path of ear {ah) between e and /. By Lemma 
2, we can replace the line segment ef (C piPi^i) by ef to obtain a new path from pi to such 
that the new path has the same Li length as piPi^. Further, as a portion of the obstacle path of 
ear {ah), ef is a boundary portion of the obstacle A' that contains ear {ah), and thus ef does not 
intersect the interior of any obstacle in T". 
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By processing each ear whose interior is intersected by PiPi+i as above, we find a new path from 
Pi to pi+i such that the path has the same Li length as PiPi+T and the path does not intersect the 
interior of any obstacle in T". 

By processing each segment PiPi+T in TVcoreis, t) as above for i = 0, 1, . . . , fc, we obtain another 
s-t path 7r(s,t) such that the Li length of 7r(s,t) is equal to that of TTcoreis^t) and 7r(,s,t) avoids all 
obstacles in V' . Below, we show that 7r(s,t) is a shortest s-t path avoiding the obstacles in V' . 

Since each core in core{V') is contained in an obstacle in V\ the length of a shortest s-t path 
avoiding core{V') cannot be longer than that of a shortest s-t path avoiding V' . Because the length 
of 7r(s, t) is equal to that of 'Kcore{s, t) and TTcore 

{s,t) is a shortest s-t path avoiding core{V'), Tr{s,t) 

is a shortest s-t path avoiding V' . 

Note that the above discussion also provides a way to construct 7r(s,i), which can be easily 
done in 0(n) time with the help of the triangulation Tn{V'). The lemma thus follows. □ 

Since each core in core{V) is contained in an obstacle in V', the corollary below follows from 
Lemma 4 immediately. 

Corollary 1 A shortest s-t path avoiding the obstacles in V' is a shortest s-t path avoiding the 
cores in coreij^'). 

3.2 Computing a Single Shortest Path 

Based on Lemma 4, our algorithm for finding a single shortest s-t path works as follows: (1) Apply 
Mitchell's algorithm [29, 30] on core{V') to find a vertex-preferred shortest s-t path avoiding the 
cores in core(V'); (2) by Lemma 4, find a shortest s-t path that avoids the obstacles in V' . The first 
step takes 0{hlogh) time and 0{h) space since the cores in core{V') have totally 0{h) vertices. 
The second step takes 0{n) time and 0{n) space. 

Theorem 1 Given a set of h pairwise disjoint convex polygonal obstacles of totally n vertices in 
the plane, we can find an Li shortest path between two points in the free space in 0{n + h\ogh) 
time and 0{n) space. 

3.3 Computing the Shortest Path Map 

In this subsection, we compute the SPM for V' . Mitchell's algorithm [29, 30] can compute an 0(n) 
size SPM with the standard query performances in 0(n log n) time and 0(n) space. 

By applying Mitchell's algorithm [29, 30] on the core set core{V'), we can compute an 0{h) 
size SPM in 0{h\ogh) time and 0{h) space, denoted by SPM{core{V'),s). With a planar point 
location data structure [13, 26], for any query point t in the free space J^, the length of a shortest 
s-t path avoiding core{V') can be reported in Oilogh) time, which is also the length of a shortest 
s-t path avoiding V' by Lemma 4. We thus have the following result. 

Theorem 2 Given a set of h pairwise disjoint convex polygonal obstacles of totally n vertices in 
the plane, in 0{n + hlogh) time and 0{n) space, we can construct a shortest path map of size 
0{h) with respect to a source point s, such that the length of an Li shortest path between s and any 
query point in the free space can be reported in Oilogh) time. 

The result in Theorem 2 is superior to Mitchell's algorithm [29, 30] in three aspects, i.e., the 
preprocessing time, the SPM size, and the length query time. However, with the SPM for Theorem 
2, an actual shortest path avoiding V' between s and a query point t cannot be reported in additional 
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time proportional to the number of turns of the path, although we can use this SPM to report an 
actual shortest path TTcorels,*) between s and t avoiding core{V') in additional time proportional 
to the number of turns of 7Tcore{s,t) and then find an actual shortest path avoiding V between s 
and t in another 0(n) time using 7Tcore{s,t) by Lemma 4. 

To process queries on actual shortest paths avoiding V' efficiently, in Lemma 5 below, using 
SPM(core{'P'),s), we compute an SPM for V' , denoted by SPM{Ai), of 0{n) size, which has 
the standard query performances, i.e., answers a shortest path length query in O(logn) time and 
reports an actual path in additional time proportional to the number of turns of the path. 

Lemma 5 Given the shortest path m,ap SPM {cor €{1^'), s) for the core set core(V'), we can compute 
a shortest path map SPM^M) for the obstacle set V in 0{n) time (with the help of the triangulation 
Tn{V')). 

Proof: Note that the polygons in V' are pairwise disjoint in their interior. For simplicity of 
discussion in this proof, we assume that any two different polygons in V' have disjoint interior as 
well as disjoint boundaries. 

Consider a cell Ccore{f) with the root r in SPM{core{V')^ s). Recall that r is always a vertex 
of a core in core{T'') and all points in Ccoreif) are visible to r with respect to core[V') [29, 30]. In 
other words, for any point p in the cell Ccore{f)i the line segment rp is contained in Ccorei^''), and 
further, there exists a shortest s-p path avoiding core{V') that contains rp. 

Denote by T{V') (resp., F{core{V'))) the free space with respect to V' (rcsp., core{V')). Note 
that the cell Ccore{f) is a simple polygon in T{core{V')). Wc assume that Ccore{f) contains some 
points in FiV') since otherwise we do not need to consider Ccore{^)- 

The cell Ccore{r) may intersect some ears. In other words, certain space in Ccoreif) may be 
occupied by some ears. Let C{r) be the subregion of Ccore{f) by removing from Ccore{f) the space 
occupied by all ears except their obstacle paths. Thus C(r) lies in T{V'). However, for each point 
p C{r), p may not be visible to r with respect to V' . Our task here is to further decompose C(r) 
into a set of SPM regions such that each such region has a root visible to all points in the region 
with respect to V'; further, we need to make sure that each point q in an SPM region has a shortest 
path in J-(V') from ,s that contains the line segment connecting q and the root of the region. For 
this, we first show that C(r) is a connected region. 

To show that C(r) is connected, it suffices to show that for any point p € C(r), there is a path 
in C(r) that connects r and p. Consider an arbitrary point p £ C{r). Since p € Ccx>re{f)^ rp is in 
Ccore{f) and there is a shortest path in T{core{V')) from s to p that contains rp. If the segment rp 
does not intersect the interior of any ear, then we are done since rp is in C(r). If rp intersects the 
interior of some ears, then let ear{ab) be one of such ears. By the proof of Lemma 4, rp penetrates 
ear{ah). Let e and / be the two points on the obstacle path of ear {ah) that rp intersects, and e/ 
be the portion of the obstacle path between e and /. Note that if rp is horizontal or vertical, then 
it cannot penetrate ear{ab) due to the monotonicity of its obstacle path by Observation 1. WLOG, 
assume rp is positive-sloped. Then by Lemma 2, ear{ab) is also positive-sloped. Recall that e and 
/ lie on rp. WLOG, assume r is higher than p and / is higher than e. Then the segment ef from e 
to / is monotonically increasing in both the x- and y-coordinates. By Observation 1, the obstacle 
path portion ef from e to / is also monotonically increasing in both the x- and y-coordinates. As 
in the proof of Lemma 4, for any point q G ef, there is a shortest path in T{core{V')) from s to 
q that contains rf and the portion of ef between / and q. Since ef is on rp contained in the cell 
Ccore{f)j by the properties of the shortest path map SPM{^A) [29, 30], ef is also contained in the 
cell Ccoreif)- Thus, ef is also contained in C{r). If we process each ear whose interior intersects rp 
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as above, we find a path in C{r) that connects r and p; further, this path has the same Li length 
as rp. Hence, C(r) is a connected region. 

Next, we claim that for any point p e C{r), there is a shortest path in T{V') from s to p that 
contains r. Indeed, since p G Ccore{r), there is a shortest path in JF{core{V')) from s io p that 
contains rp; let T^coreis-, r) be the portion of this path between s and r. On one hand, wc have shown 
above that there is a path from r to p in C(r) with the same Li length as rp. On the other hand, 
by Lemma 4, there exists a path in F{V') from s to r with the same length as -Kcoreis^r). Hence, 
a concatenation of these two paths results in a shortest path from s to p in F{V') that contains r. 
Our claim thus follows. 

The above claim and its proof also imply that decomposing C(r) into a set of SPM regions is 
equivalent to computing an SPM in C(r) with the vertex r as the source point, which we denote by 
SPM{C{r)). Since C(r) is a connected region and Ccorei^^ is a simple polygon, we claim that C(r) 
is a (possibly degenerate) simple polygon. This is because for any ear E that intersects Ccore{r), 
the portion E H Ccoreif) lies on the boundary of the simple polygon Ccoreif)', thus, removing E 
except its obstacle path from Ccoreif) (to form C(r)) changes only the boundary shape of Ccore{f) 
but does not change the nature of a simple polygonal region (from Ccore{i') to C(r)). Based on 
the fact that C(r) is a (possibly degenerate) simple polygon, SPM{C{r)) can be easily computed 
in linear time in terms of the number of edges of C(r). For example, since the Euclidean shortest 
path between any two points in a simple polygon is also an Li shortest path between the two points 
[17], an SPM in a simple polygon with respect to the Euclidean distance is also one with respect to 
the Li distance. Therefore, we can use a corresponding shortest path algorithm for the Euclidean 
case (e.g., [16]) to compute each SPM{C{r)) in our problem. 

Note that our discussion above also implies that given SPM{core{V')^ s), for each cell Ccore{f) 
with a root r, we can compute the corresponding SPM{C{r)) separately. Clearly, the SPM{C{r))'s 
corresponding to all cells in SPM {core{V' ) , s) constitute a shortest path map SPM{M) for V'. 

Due to the planarity of the cell regions involved, the total number of edges of all C(r)'s is 0(n). 
Given a triangulation Tn{'P'), all regions C(r) can be obtained in totally 0(n) time. Computing 
all SPM{C{r)ys also takes totally 0(n) time. Thus, SPM{M) can be constructed in 0{n) time. 
The lemma thus follows. □ 

Theorem 2 and Lemma 5 together lead to the following result. 

Theorem 3 Given a set of h pairwise disjoint convex polygonal obstacles of totally n vertices in 
the plane, in 0{n + hlogh) time and 0{n) space, we can construct a shortest path map of size 0{n) 
with respect to a source point s, such that given any query point t in the free space, the length of an 
Li shortest s-t path can be reported in 0(logh) time and an actual path can be found in 0{logn+k) 
time where k is the number of turns of the path. 

4 Shortest Paths among General Polygonal Obstacles 

In this section, we consider the general case, i.e., the obstacles in V are not necessarily convex. In 

the following, in Section 4.1, we review the corridor structure [25], and introduce the ocean A4. In 
Section 4.2, we present the algorithm for computing a single shortest path and the similar idea also 
computes an SPM for M, i.e., SPM{^A). In Section 4.3, we outline our algorithm for computing 
an SPM in the entire free space T. 
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Figure 4: Illustrating a triangulation of the free space among 
two obstacles and the corridors (with red solid curves). There 
are two junction triangles indicated by the large dots inside 
them, connected by three solid (red) curves. Removing the 
two junction triangles results in three corridors. 



Figure 5: Illustrating an open hourglass (left) and 
a closed hourglass (right) with a corridor path linking 
the apices x and y of the two funnels. The dashed 
segments are diagonals. The paths 7r(a, b) and n{e, /) 
are shown with thick solid curves. 



4.1 Preliminaries 

For simplicity of discussion, we assume that all obstacles are contained in a large rectangle TZ (see 
Fig. 4). Let J- be the free space inside TZ. Let t be an arbitrary point in J^. 

We first review the corridor structure [25]. Denote by Tn{J-) a triangulation of J-". Let G{J-) 
denote the (planar) dual graph of Tn{J^), i.e., each node of G{J^) corresponds to a triangle in Tn{T) 
and each edge connects two nodes of G{J-) corresponding to two triangles sharing a diagonal of 
Tn{T). The degree of each node in G{F) is at most three. As in [25], at least one node dual to a 
triangle incident to each of s and t is of degree three. Based on G{F), we compute a planar 3-regular 
graph, denoted by G^ (the degree of each node in G^ is three), possibly with loops and multi-edges, 
as follows. First, we remove every degree-one node from G{J-) along with its incident edge; repeat 
this process until no degree-one node exists. Second, remove every degree-two node from G{J^) 
and replace its two incident edges by a single edge; repeat this process until no degree-two node 
exists. The resulting graph is G^ (e.g., see Fig. 4). The resulting graph G^ has 0{h) faces, 0{h) 
nodes, and 0{h) edges [25]. Each node of G^ corresponds to a triangle in Tri{J-), which is called 
a junction triangle (e.g., see Fig. 4). The removal of all junction triangles from G^ results in 0{h) 
corridors, each of which corresponds to one edge of G^. 

The boundary of a corridor G consists of four parts (see Fig. 5): (1) A boundary portion of 
an obstacle Pi G V, from a point a to a point b; (2) a diagonal of a junction triangle from 6 to a 
boundary point e on an obstacle Pj G V {Pi = Pj is possible); (3) a boundary portion of the obstacle 
Pj from e to a point /; (4) a diagonal of a junction triangle from / to a. The two diagonals be and af 
are called the doors of C. The corridor C is a simple polygon. Let 7r(a, b) (resp., 7r(e, /)) denote the 
shortest path from a to 6 (resp., e to /) inside C. The region He bounded by 7r(a, 6), 7r(e, /), and 
the two diagonals be and fa is called an hourglass, which is open if 7r(a, b) D 7r(e, f) = $ and closed 
otherwise (see Fig. 5). If He is open, then both 7r(a, b) and 7r(e, /) are convex chains and are called 
the sides of He; otherwise, He consists of two "funnels" and a path ttc = vr(o, b) fl 7r(e, /) joining 
the two apices of the two funnels, called the corridor path of G. The two funnel apices connected 
by the corridor path are called the corridor path terminals. Each funnel side is also convex. We 
compute the hourglass for each corridor. After the triangulation, computing the hourglasses for all 
corridors takes totally 0[n) time. 

Let Q be the union of all junction triangles and hourglasses. Then Q consists of 0{h) junction 
triangles, open hourglasses, funnels, and corridor paths. As shown in [21], there exists a shortest 
s-t path 7r(s, t) avoiding the obstacles in V which is contained in Q. Consider a corridor C. If 
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7r(s, t) contains an interior point of C, then the path 7r(s, t) must intersect both doors of C; further, 
if the hourglass Hq of C is closed, then we claim that we can make the corridor path of C entirely 
contained in 7r(s,t). Suppose 7r(s,t) intersects the two doors of C, say, at two points p and q 
respectively. Then since C is a simple polygon, a Euclidean shortest path between p and q inside 
C, denoted by 7Te{p,q), is also an Li shortest path in C [17]. Note that 7Te{p,q) must contain the 
corridor path of C. If we replace the portion of 7r(s, t) between p and q by tte{p^ q), then we obtain 
a new Li shortest s-t path that contains the corridor path ttc- For simplicity, we still use 7r{s,t) to 
denote the new path. In other words, Tr{s,t) has the property that if 7r(s,t) intersects both doors 
of C and the hourglass He is closed, then the corridor path of C is contained in 7r(s,t). 

Let A4 he Q minus the corridor paths. We call A4 the ocean. Clearly, A4 Q T. The boundary 
of Ai consists of 0{h) reflex vertices and 0{h) convex chains, implying that the complementary 
region Tl\M consists of a set of polygons of totally 0{h) reflex vertices and 0{h) convex chains. 
As shown in [25], the region Tl\M. can be partitioned into a set V of 0{h) convex polygons of 
totally 0{n) vertices (e.g., by extending an angle-bisecting segment inward from each reflex vertex). 
The ocean A4 is exactly the free space with respect to the convex polygons in V' . In addition, for 
each corridor path, no portion of it lies in M. Further, the shortest path iT{s,t) is a shortest s-t 
path avoiding all convex polygons in V' and possibly utilizing some corridor paths. The set T" can 
be easily obtained in 0(n + hlogh) time. Therefore, as in [25], other than the corridor paths, we 
reduce our original Li-SP problem to the convex case. 

4.2 Finding a Single Shortest Path and Computing an SPM for M 

With the convex polygon set V' , to find a shortest s-t path in J^, if there is no corridor path, 
then we can simply apply our algorithm for the convex case in Section 3. Otherwise, the situation 
is more complicated because the corridor paths can give possible "shortcuts" for the sought s-t 
path, and we must take these possible "shortcuts" into consideration while running the continuous 
Dijkstra paradigm [29, 30]. The details are given below. 

First, we compute the core set core{V') of V' . However, the way we construct core(V') here 
is slightly different from Section 3. For each convex polygon A' e T", in addition to its leftmost, 
topmost, rightmost, and bottommost vertices, if a vertex v of A' is a corridor path terminal, then v 
is also kept as a vertex of the core core{A'). In other words, core{A') is a simple (convex) polygon 
whose vertex set consists of the leftmost, topmost, rightmost, and bottommost vertices of A' and 
all corridor path terminals on A'. Since there are 0{h) terminal vertices, the cores in core{V') still 
have totally 0{h) vertices and edges. Further, the core set thus defined still has the properties 
discussed in Section 3 for computing shortest Li paths, e.g., Observation 1 and Lemmas 1, 2, and 
4. Hence, by using our scheme in Section 3, we can first find a shortest s-t path avoiding the cores 
in core{V') in O {hlogh) time by applying Mitchell's algorithm [29, 30], and then obtain a shortest 
s-t path avoiding V' in 0{n) time by Lemma 4. But, the path thus computed may not be a true 
shortest path in J- since the corridor paths are not utilized. To find a true shortest path in we 
need to modify the continuous Dijkstra paradigm when applying it to core{V'), as follows. 

In Mitchell's algorithm [29, 30], when an obstacle vertex v is hit by the wavefront for the first 
time, it will be "permanently labeled" with a value d{v), which is the length of a shortest path from 
s to f in the free space. The wavefront consists of many "wavelets" (each wavelet is a line segment of 
slope 1 or —1). The algorithm maintains a priority queue (called "event queue"), and each element 
in the queue is a wavelet associated with an "event point" and an "event distance" , which means 
that the wavelet will hit the event point at the event distance. The algorithm repeatedly takes 
(and removes) an element from the event queue with the smallest event distance, and processes 
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the event. After an event is processed, some new events may be added to the event queue. The 
algorithm stops when the point t is hit by the wavefront for the first time. 

To handle the corridor paths in our problem, consider a corridor path ttq with x and y as its 
terminals and let I be the length of ttc- Recall that x and y are vertices of a core in core{V'). 
Consider the moment when the vertex x is permanently labeled with the distance d{x). Suppose 
the wavefront that first hits x is from the funnel whose apex is x. Then according to our discussions 
above, the only way that the wavelet of the wavefront at x can affect a shortest s-t path is through 
the corridor path ttc- If y is not yet permanently labeled, then y has not been hit by the wavefront. 
We initiate a "pseudo-wavelet" that originates from x with the event point y and event distance 
d{x) + /, meaning that y will be hit by this pscudo- wavelet at the distance d{x) + /. We add the 
pseudo-wavelet to the event queue. If y has been permanently labeled, then the wavefront has 
already hit y and is currently moving along the corridor path ttc from y to x. Thus, the wavelet 
through X will meet the wavelet through y somewhere on the path ttc, and these two wavelets will 
"die" there and never affect the free space outside the corridor. Thus, if y has been permanently 
labeled, then we do not need to do anything on y. In addition, at the moment when the vertex x 
is permanently labeled, if the wavefront that first hits x is from the corridor path ttc (i-e., through 
y), then the wavelet at x will keep going to the funnel of x through x; therefore, we process this 
event on x as usual (i.e., as in [29, 30]), by initiating new wavelets that originate from x. 

For a corridor path ttc with two terminals x and y, when x is permanently labeled, if the 
wavefront that first hits x is not from the corridor path ttc, then we call x a wavefront incoming 
terminal; otherwise, a; is a wavefront outgoing terminal. According to our discussion above, at least 
one of X and y must be a wavefront incoming terminal. In fact, both x and y can be wavefront 
incoming terminals, in which case the wavefronts passing through x and y "die" inside the corridor. 

Intuitively, the above treatment of corridor path terminals makes corridor paths act as possible 
"shortcuts" when we propagate the wavefront. The rest of the algorithm proceeds in the same 
way as in [29, 30] (e.g., processing the segment dragging queries). The algorithm stops when the 
wavefront first hits the point t, at which moment a shortest s-t path in T has been found. 

Since there are 0{h) corridor paths, with the above modifications to Mitchell's algorithm as 
applied to coreiV'), its running time is still 0{h\ogh). Indeed, comparing with the original con- 
tinuous Dijkstra scheme [29, 30] (as applied to core{P')), there are 0{h) additional events on the 
corridor path terminals, i.e., events corresponding to those pseudo- wavelets. To handle these addi- 
tional events, we may, for example, as preprocessing, for each corridor path, associate with each its 
corridor path terminal x the other terminal y as well as the corridor path length I. Thus, during 
the algorithm, when we process the event point at x, we can find y and I immediately. In this way, 
each additional event is handled in 0(1) time in addition to adding a new event for it to the event 
queue. Hence, processing all events still takes 0{h\ogh) time. Note that the shortest s-t path thus 
computed may penetrate some ears of V' . As in Lemma 4, we can obtain a shortest s-t path in the 
free space T in additional 0(n) time. Since applying Mitchell's algorithm on coreiV') takes 0{h) 
space, the space used in our entire algorithm is 0(n). 

In summary, we have the following result. 

Theorem 4 Given a set of h pairwise disjoint polygonal obstacles of totally n vertices in the plane, 
we can find an Li shortest path between two points in the free space in 0{n -\- /ilog^^^ h) time (or 
0{n-\- h\og h) time if a triangulation of the free space is given) and 0{n) space. 

As Mitchell's algorithm [29, 30], the above algorithm also computes a shortest path map on the 
free space of the convex polygons in V, i.e., SPM{Ai). We should point out that because of the 
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Figure 6: Illustrating a bay bay{cd) in an open hourglass (left) and a canal canal{x, y) in a closed hourglass (right) 
with a corridor path linking the apices x and y of its two funnels. 

0{h) corridor paths, SPM{A4) is different from a "normal" SPM in the fohowing aspect. Consider 
a corridor path ttc with two terminals x and y. Suppose x is a wavefront incoming terminal and 
y is a wavefront outgoing terminal. Then this means that the algorithm determines a shortest 
path from s to y which goes through x. Corresponding to the corridor path ttc, we may put a 
"pseudo-cell" in SPM{A4) with x as the root such that y is the only point in this "pseudo-cell", 
and we also associate with the pseudo-cell the corridor path ttc, which indicates that there is a 
shortest s-y path that consists of a shortest s-x path and the corridor path nc- If x and y are 
both wavefront incoming terminals, then we need not do anything for this corridor path. Clearly, 
since there are 0{h) corridor paths, the above procedure of building pseudo-cells affects neither the 
space bound nor the time bound for constructing SPM{Ai). Therefore, the SPM{A4) of size 0(n) 
can be computed in 0{T) time and 0{n) space, where T is the time for triangulating T. Based on 
SPM{A4), in Section 4.3, we will compute an SPM on the entire free space in additional 0{n) 
time. 

4.3 Computing a Shortest Path Map 

Based on SPM{A4), in Section 4.3, together with Sections 5 and 6, we will compute in additional 
0{n) time an SPM on the entire free space with respect to the source point s, denoted by 
SPM{J^), which has the standard query performances, i.e., for any query point t, it reports the 
length of a shortest s-t path in O(logn) time and the actual path in additional time proportional 
to the number of turns of the path. 

As discussed in [29, 30], SPM{T) may not be unique. We show that an SPM{F) of size 0{n) 
can be computed in 0{n + /ilog^^*^ h) time (or 0{n -\- hlogh) time if a triangulation of the free 
space is given). Our techniques for constructing SPM {J-") are independent of those in the earlier 
sections of this paper, and are also different from those in the previous work (e.g., [29, 30]). 

This section introduces the new concepts, bays and canals, and outlines the algorithm, while 
the details are given in Sections 5 and 6. One key subproblem we need to solve efficiently is the 
special weighted Li geodesic Voronoi diagram problem, i.e., the challenging subproblem illustrated 
in Fig. 1. A linear time algorithm is given in Section 5 for it. Section 6 deals with another 
subproblem, where the algorithm in Section 5 is used as a procedure. 

4.3.1 Bays and Canals 

Recall that Ai T. To compute SPM{F), since we already have SPM{M.), we only need to 
compute the portion of SPM{J^) in the space J-\M.. We first examine the space J-\M, which we 
partition into two type of regions, hays and canals, defined as follows. 
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Consider an hourglass He of a corridor C. We first discuss the case when He is open (see Fig. 6). 
He has two sides. Let Si{Hc) be an arbitrary side of He- The obstacle vertices on Si{Hc) all lie 
on the same obstacle, say P G "P. Let c and d be any two adjacent vertices on Si{Hc) such that 
the line segment cd is not an edge of P (see the left figure in Fig. 6, with P = Pj). The region 
enclosed by cd and a boundary portion of P between c and d is called the bay of P and cd, denoted 
by bay{cd), which is a simple polygon. We call cd the bay gate. 

If the hourglass He is closed, then let x and y be the two apices of its two funnels. Consider 
two adjacent vertices c and d on a side of a funnel such that the line segment cd is not an obstacle 
edge. If neither c nor d is a funnel apex, then c and d must both lie on the same obstacle and the 
segment cd also defines a bay with that obstacle as above. However, if either c or d is a funnel 
apex, say, x = c, then x and d may lie on different obstacles. If they both lie on the same obstacle, 
then they also define a bay; otherwise, we call xd the canal gate at x (see Fig. 6). Similarly, there 
is also a canal gate at the funnel apex y, say yz. Let Pi and Pj be the two obstacles defining the 
hourglass He- The region enclosed by Pi, Pj, and the two canal gates xd and yz that contains the 
corridor path of He is called the canal of He, denoted by canal{x,y), which is a simple polygon. 

It is easy to see that J^\M. consists of all bays and canals thus defined. 

To build SPM{T), we need to compute the portion of SPM{F) in all bays and canals since 
we already have SPM(Ai). As all bays and canals are connected with A4 through their gates, 
we need to "expand" SPM{Ai) to all bays/canals through their gates. Henceforth, when saying 
"compute an SPM for a bay/canal," we mean "expand SPM{AA) into that bay/canal", and vice 
versa. Computing an SPM for a bay is a key (i.e., the challenging subproblem). Computing an 
SPM for a canal uses the algorithm for a bay as a main procedure. 

4.3.2 Expanding SPM{M) into Bays and Canals 

We discuss the bays first. Consider a bay bay{cd). If its gate cd is in a single cell C(r) of SPM{A4) 
with r as the root, then each point in bay{cd) has a shortest path to s via r. Thus, to construct 
an SPM for bay{cd), it suffices to compute an SPM on bay{cd) with respect to the single point r. 
This can be easily done in linear time (in terms of the number of vertices of bay{cd)) since bay{cd) 
is a simple polygon*. Note that although r may not be a vertex of hay{cd), we can, for example, 
connect r to both c and d with two line segments (both fc and rd are in C{r)) to obtain a new 
simple polygon that contains bay{cd). 

If the gate cd is not contained in a single cell of SPM{M), then the situation is more complicated. 
In this case, multiple vertices of SPM(A4) may lie in the interior of cd (i.e., the intersections of 
the boundaries of the cells of SPM{M) with cd). This is actually the challenging subproblem 
illustrated by Fig. 1. We refer to the vertices of SPM{M.) on cd (including its endpoints c and d) 
as the SPM{M.) vertices and let m! be their total number. Let n! be the number of vertices of 
bay{cd). A straightforward approach for computing an SPM for bay{cd) is to use the continuous 
Dijkstra paradigm [29, 30] to let the wavefront continue to move into bay{cd). But, this approach 
may take 0{{n' +m') log(m' + n')) time. Later in Section 5, we derive an 0{n' + m') time algorithm, 
as stated below. 

Theorem 5 For a bay of n' vertices with m' SPM{M) vertices on its gate, a shortest path map of 
size 0{n' + m!) for the bay can be computed in 0{n' + m') time. 

*For example, since the Euclidean shortest path between any two points in a simple polygon is also an Li shortest 
path [17], a Euclidean SPM in a simple polygon is also an Li one. Thus, we can use a corresponding shortest path 
algorithm for the Euclidean case (e.g., [16]) to compute an Li SPM in bay{cd) with respect to r in linear time. 
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Since a canal has two gates which are also edges of M., multiple SPM{^A) vertices may lie on 
both its gates. Later in Section 6, we show the following result. 

Theorem 6 For a canal of n' vertices with totally m! SPM{M) vertices on its two gates, a shortest 
path map of size 0{n' + m') for the canal can be computed in 0{n' + m') time. 

4.3.3 Wrapping Things Up 

By Theorems 5 and 6, the time bound for computing the shortest path maps for all bays and canals 
is linear in terms of the total sum of the numbers of obstacle vertices of all bays and canals, which 
is 0(n), and the total number of the SPM{A4) vertices on the gates of all bays and canals, which 
is also 0{n) since the size of SPM{M.) is 0{n). 

We hence conclude that given SPM {A4), SPM{T) can be computed in additional 0(n) time. 
With a linear size planar point location data structure [13, 26], we have the following result. 

Theorem 7 Given a set of h pairwise disjoint polygonal obstacles of totally n vertices and a 
source point s in the plane, we can build a shortest path map of size 0{n) with respect to s in 
0(n + h\og^~^'^ h) time (or 0(n + h\ogh) time if a triangulation of the free space is given) andO{n) 
space, such that for any query point t, the length of a shortest s-t path can be reported in O(logn) 
time and the actual path can be found in additional 0{k) time, where k is the number of turns of 
the path. 

5 Computing a Shortest Path Map for a Bay 

Consider a bay hoy(cd) with the gate cd (see Fig. 6). Let SPM (bay (cd)) be the SPM for bay{cd) 
that we seek to compute. 

For the case when the segment cd lies in a single cell C(r) of SPM{M.) with the root r, we have 
already shown how to construct SPM{bay{cd)) in linear time (in terms of the number of vertices 
of bay{cd)). If the gate cd is not contained in a single cell of SPM{A4), then let m' be the number 
of SPM(Ai) vertices on cd, and n' be the number of vertices of bay{cd). In this section, we give an 
algorithm for computing SPM{bay{cd)) in 0{n' + m!) time. 

Let R be the set of roots of the cells of SPM{M.) that intersect with cd. To obtain SPM{bay{cd)), 
we can first compute, for each r ^ R, the Voronoi region VD{r) inside bay{cd) such that for any 
point t G VD{r), there is a shortest s-t path via r; we then compute an SPM on VD{r) with respect 
to the single point r. Since every VD{r) is a simple polygonal region in bay{cd), the shortest path 
map SPM{VD{r),r) can be computed in linear time in terms of the number of vertices of VD{r) 
(e.g., by using an algorithm in [16, 17]). Thus, the key is to decompose bay{cd) into Voronoi regions 
for the roots of R, which is exactly the challenging subproblem illustrated by Fig. 1. Denote by 
VD{bay{cd)) this Voronoi diagram decomposition of bay{cd). We aim to compute VD{bay{cd)) in 
0{n' + m') time. 

Without loss of generality (WLOG), assume that cd is positive-sloped, bay{cd) is on the right 
of cd, and the vertex c is higher than d (e.g., bay{cd) = B va. Fig. 1). Other cases can be handled 
similarly. Let R = {ri,r2, . . . ,rk} be the set of roots of the cells of SPM[JV[) that intersect with 
cd in the order from c to d along cd. Note that R may be a multi-set, i.e., two roots rj and Vj with 
i ^ j may refer to the same physical point; but this is not important to our algorithm (e.g., we 
can view each as a physical copy of the same root). Let c = vq,vi, . . . ,Vk = d be the SPM{M.) 
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vertices on cd ordered from c to d (thus m' = A; + 1). Hence, for each 1 < i < k, the segment 
Vi-iVi is on the boundary of the ceh C(rj) of SPM(M). Note that each cell C(rj) is a star-shaped 
polygon, and for each I < i < k — 1, Vi lies on the common boundary of C{ri) and C(ri+i) (i.e., 
Vi G C{ri) n C(ri+i)). To obtain VD{bay{cd)), for each Vi G R, we need to compute the Voronoi 
region VD(rj). 

Our algorithm can be viewed as an incremental one, i.e., it considers the roots in R one by one. 
It is commonly known that incremental approaches can construct Voronoi diagrams in quadratic 
time, or may give good randomized result. In contrast, our algorithm is deterministic and takes 
only linear time. The success of it hinges on that we can find an order of the roots in R such that by 
following this order to consider the roots in R incrementally, we are able to compute VD{bay{cd)) 
in linear time. The order is nothing but that of the indices of the roots in R we have defined. With 
this order, the algorithm is quite simple. However, it is quite challenging to argue its correctness 
and achieve a linear time implementation. Our strategy is to show that the algorithm implicitly 
maintains a number of invariants that assure the correctness of the algorithm. For this purpose, 
we give many observations (in Section 5.2). Additionally, some interesting techniques are also used 
to implement and simplify the algorithm. 

We first give an algorithm overview in Section 5.1. 

5.1 Algorithm Sketch 

To compute VD(bay{cd)), it turns out that we need to deal with the interactions between some 
rays, each of which belongs to the bisector of two roots in R. Every such ray is either horizontal 
or vertical. Further, considering the roots in R incrementally is equivalent to considering the 
corresponding rays incrementally. We process these rays in a certain order (e.g., as to be proved, 
their origins somehow form a staircase structure). For each ray considered, if it is vertical, then it 
is easy (it eventually leads to a ray shooting operation), and its processing does not introduce any 
new ray. But, if it is horizontal, then the situation is more complicated since its processing may 
introduce many new horizontal rays and (at most) one vertical ray, also in a certain order along 
a staircase structure (in addition to causing a ray shooting operation). A stack is used to store 
certain vertical rays that need to be further processed. 

The algorithm needs to perform ray shooting operations for some vertical and horizontal rays. 
Although there are known data structures for ray shooting queries [4, 5, 16, 18], they are not 
efficient enough for a linear time implementation of the entire algorithm. Based on observations, 
our approach makes use of the horizontal visibility map and vertical visibility map of bay{cd) [3]. 
More specifically, we prove that all vertical ray shootings are in a "nice" sorted order (called target- 
sorted). With this property, all vertical ray shootings are performed in totally linear time by using 
the vertical visibility map of hay{cd). The horizontal visibility map is used to guide the overall 
process of the algorithm. During the algorithm, we march into the bay and the horizontal visibility 
map allows us to keep track of our current position (i.e., in a trapezoid of the map that contains 
our current position). The horizontal visibility map also allows each horizontal ray shooting to be 
done in 0(1) time. In addition, in the preprocessing of the algorithm, we also need to perform 
some other ray shootings (for rays of slope —1); our linear time solution for this also hinges on the 
target-sorted property of such rays. 

Our algorithm is conceptually simple. As mentioned above, the only data structures we need 
are linked lists, a stack, and the horizontal and vertical visibility maps. Its correctness relies on 
the fact that the algorithm implicitly maintains a set of invariant properties in each iteration. To 
prove the algorithmic correctness, of course, we need to show that these invariant properties hold 
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Figure 7: Illustrating some cases of the bisector -B(pi,p2) of two weighted points pi and p2. In (3), an entire 
quadrant (the shaded area) can be used as B{pi,p2), but we choose B{pi,p2) to be the vertical (solid thick) half-line. 

itcrativcly. Specifically, in our discussion of the algorithm, after each iteration we formally prove 
that the invariants are well maintained. For this purpose, before presenting the algorithm in Section 
5.3, we first show a set of observations in Section 5.2, which capture some essential properties of 
this Li problem. These observations may be helpful for solving other related problems as well. 
However, the discussion of these observations and the formal proofs that the invariant properties 
are maintained by the algorithm somehow make the presentation of this whole section lengthy, 
technically complicated, or even tedious, for which we ask for the reader's patience. 



5.2 Observations 

In this subsection, we give a number of observations, most of which help capture the behaviors 
of the bisectors for the roots of R in computing VD{bay{cd)). Although some of the observations 
individually might appear simple, they are essential and adding them up leads to an efficient 
algorithmic strategy for computing VD{bay{cd)) (as presented in Section 5.3). The observations 
also allow our algorithm to perform some key operations (e.g., ray shootings) in a faster manner 
than using a standard approach [4, 5, 16, 18]. 

For a point p, denote by x{p) its x-coordinate and by y{p) its y-coordinate. For two objects Oi 
and O2 in the plane, if x{pi) < x{p2) for any two points pi G 0\ and p2 € O2, then we say Oi is to 
the left or west of O2, or O2 is to the right or east of Oi\ if y{pi) < y{p2) for any two points pi € Oi 
and P2 & O2, then we say Oi is to the south of O2 or Oi is below O2 , or O2 is to the north of Oi or 
O2 is above 0\. If Oi is to the left of O2 and is also below O2, then we say 0\ is to the southwest 
of O2 or O2 is to the northeast of Oi. We define southeast and northwest similarly. 

In our problem, each root ri ^ R can be viewed as an additively weighted point whose weight 
is the Li length of a shortest path from s to r^. Thus, we need to consider the possible shapes 
of the bisector of two weighted points. For two weighted points pi and p2 with weights wi and 
W2-, respectively, their bisector B{pi,p2) consists of all points q such that the Li length of the line 
segment piq plus wi is equal to the Li length of p2q plus W2- Figure 7 shows some cases. Note 
that the bisector can be an entire quadrant of the plane (e.g., see Figure 7(3)); in this case, as in 
[29, 30], we choose a vertical half-line as the bisector. For any pair of consecutive roots r^-i and 
rj in i? for 2 < i < k, since the SPM{M) vertex Vi-i G cd is on the common boundary of C(ri_i) 
and C{ri), Vi-i lies on the bisector B{ri^i,ri) of r.j_i and r^. For two points pi and p2, denote 
by Rec{pi,p2) the rectangle with pi and p2 as its two diagonal vertices. The next observation is 
self-evident. 

Observation 2 The bisector B{pi,p2) consists of three portions: Two half-lines and a line segment 
connecting them; the line segment has a slope 1 or —1 and is the intersection of B{pi,p2) and the 
rectangle Rec{pi,p2), and each of the two half-lines is perpendicular to an edge of Rec{pi,p2) that 
touches the half-line. Depending on the relative positions and weights of pi and p2, some portions 
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Figure 8: An example of rt to the northeast of Vi-i. The point p / Vi-i is on Vi-iVi and is infinitely close to Vi-i. 
The line segment Vip must cross n-iVi-i. 

of B{pi,p2) may degenerate and become empty. B{p\,p2) is monotone to both the x- and y-axes. 
For any line I containing a portion of B{pi,p-2), pi and p2 cannot lie strictly on the same side of I. 

We call the open line segment of B{pi,p2) strictly inside Rec{pi,p2) its middle segment, denoted 
by Bm{pi,P2), and the two half- lines of B{pi,p2) its two rays, each originating at a point on an edge 
of Rec{pi.p2). Thus, the origins of the two rays of B{p\,p2) arc the two cndpoints of B\i{pi,p2). 

Since each cell in an SPM is a star-shaped simple polygon, the observation below is obvious. 

Observation 3 Let C{r) and C{r') be two different cells in SPM{M) with roots r and r' . For any 
two points p E C{r) and p' G C{r'), the line segments pf andp'r' cannot cross each other. 

The next lemma shows the possible relative positions of two consecutive roots in R. 

Lemma 6 For any two consecutive roots rj_i and ri in R with 2 < i < k, ri cannot be to the 
northeast of ri-i, or equivalently, ri-\ cannot be to the southwest of ri. 

Proof: Since the SPM{M.) vertex Vi-i G cd lies on the common boundary of the two cells C(rj_i) 

and C{ri), Vi^i is on the bisector B{ri^i,ri). 

Assume to the contrary that is to the northeast of rj-i. Note that Vi-i may lie on either a 
half- line or the middle segment of B{ri-i,ri). In either case, since is to the northeast of rj-i and 
cd is positive-sloped, according to Observation 2, vi-i must be lower than rj, and Vj-i must be to 
the right of rj-i (see Fig. 8). 

Since the segment Vi^iVi is not a single point and vi is to the left of Vi^i, we can find a point 
p G Vi-iVi such that p ^ Vi-i and p is infinitely close to Vi-i (see Fig. 8). Since p G Vi-iVi and 
Vi-iVi C C{ri), we have p G C(rj). Note that Vi-i G C(ri_i) fl C{ri). Below we show that the two 
line segments r^p and rj_iUj_i must cross each other, which contradicts with Observation 3. 

Since both rj and rj_i are obstacle vertices, by our assumption, r.j and rj_i do not lie on a 
horizontal or vertical line. Hence is strictly to the northeast of ri_i. Note that no root in R 
lies on cd. Since Vi-i is lower than and is to the right of r^-i, the three points Vi-i, ri, and 
Tj-i do not lie on the same line (see Fig. 8). In other words, the triangle AriVj-irj-i is a proper 
one. Further, suppose p{ri,Vi-i) (resp., p(rj,rj_i)) is the ray originating from and going through 
Vi-i (resp., Tj-i); then p(ri,ri^i) can be obtained by rotating p{ri,Vi-i) clockwise by an angle 
Zvi-iriri-i > 0°. By the definition of the point p, during this rotation, p will be encountered by 
the rotating ray p{ri,Vi-i) at an angle Zvi-irip with 0° < Zvi-irip < Zvi-iriri-i, which implies 
that Tip crosses rj_if The lemma thus follows. □ 

By Lemma 6, there arc three cases on the possible relative positions of r^-i with respect to r,, 
i.e., Tj-i can be to the southeast, northwest, or northeast of r^. 

Lemma 7 Consider any two consecutive roots r^-i and ri in R with 2 < i < k. 
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Figure 9: Illustrating the three possible relative positions of ri-i and rj. 

1. If ri is to the southeast of ri-\, then ■Uj-i is on a ray of B{ri-i,ri) that is horizontally going 
east and Vi-i is to the right of Rec{ri-i,ri) (see Fig. 9(1)). 

2. If ri is to the northwest of ri-i, then Vi-i is on a ray of B{ri-i,ri) that is vertically going 
south and Vi-i is below Rec{ri-i,ri) (see Fig. 9(2)). 

3. If ri is to the southwest of ri-\, then Vi-\ is either on the middle segment BM{ri-i,ri), or 
on a ray of B{ri-i,ri) that is either horizontally going east or vertically going south (see 
Fig. 9(3)). Further, i/f j-i is on the ray horizontally going east, then Vi-i is to the right of 
Rec{ri-i,ri); if Vi-i is on the ray vertically going south, then Vi-i is below Rec{ri-\,ri). 

Proof: We first prove Part 1 of the lemma. If rj is to the southeast of r^-i (see Fig. 9(1)), then 
the rectangle i?ec(rj_i, r^) cannot intersect cd. Thus, Vi^i cannot be on BM{ri-i,ri), and Vi^i 
must be on a ray of B(ri-i,ri), denoted by p. By Observation 2, the origin of p is on an edge a of 
Rec{ri-i,ri) and is perpendicular to the edge a. Since Vi-i G p and rj-i is to the northwest of rj, a 
must be one of the two edges incident to ri, i.e., the bottom edge or the right edge of Rec{ri-i,ri). 
In addition, if a is the bottom edge of Rec{ri-i,ri), then p must be vertically going south; further, 
since is to the southeast of r-j-i, by a similar argument as that for the proof of Lemma 6, we 
can obtain a contradiction. Thus, a is the right edge of Rec{ri-i,ri) and p must be horizontally 
going east. In addition, it is easy to see that Vi-i must be to the right of Rec{ri-i,ri). Part 1 of 
the lemma thus follows. 

Part 2 can be proved analogously as Part 1, and we omit it. 

For Part 3, if Rec{ri-i,ri) intersects cd, then it is possible that BM{ri-i,ri) intersects cd (at 
Vi-\). If BM{ri~i,ri) doest not intersect cd, then Vi-i lies on a ray of B{ri-\,ri), denoted by p. 
Again, the origin of p is on either the right edge of Rec{ri-i,ri) or the bottom edge of Rec{ri-i,ri). 
In the former case, p is horizontally going east and Vi-i is to the right of iiec(rj_i, rj). In the latter 
case, p is vertically going south and Vi-i is below Rec{ri-i,ri). Part 3 thus follows. □ 

For any two consecutive roots rj_i and in R with 2 <i < k,\i Vi-i is on a ray p of B{ri-i,ri), 
then we let pi-i be the ray originating at Vi-\ with the same direction as p. If vi-i lies on the 
middle segment of B{ri-i,ri), then by Lemma 7, ri^i is to the northeast of and cd intersects 
Rec{ri-i,ri); in this case, let pi-i be the ray of B{ri^i,ri) that is below or to the right of and 
goes inside bay{cd). For a ray p, let or{p) denote the origin of p. Observation 4 below is obvious. 

Observation 4 For any 2 < i < k, the ray is either horizontally going east or vertically going 
south. If Vi-i is on a ray of B{ri-i,ri), then or{pi-i) = Vi-i; if Vi-i is on BM{ri-i,ri), then 
or{pi-i) is on either the right edge or the bottom edge of Rec{ri-i,ri). 

Lemma 8 Consider any two consecutive roots r^-i and ri in R with 2 < i < k. 
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1. If the ray pi-\ is horizontal, then rj-i is above pi-i and ri is below Pi-i- 

2. If pi-i is vertical, then rj-i is to the right of pi-i and r^ is to the left of pi-i. 

3. The origin or{pi-i) of pi-i is always below ri_i and to the right ofri. 

Proof: There axe three cases on the possible relative positions of rj_i and r^. 

• If Tj-i is to the northwest of (see Fig. 9(1)), then by the proof of Lemma 7, pi-\ is horizontal 
and is contained in the ray of B{ri-i,ri) whose origin is on the right edge of Rec{ri-i,ri). 
Since rj_i and are two diagonal vertices of Rec{ri-i,ri), pi-i is above and below rj_i. 

Further, the origin or{pi-\) is vi-i, which is below rj_i and to the right of r^. 

• If r^-i is to the southeast of r^ (see Fig. 9(2)), then by the proof of Lemma 7, pi-\ is vertical 
and lies on the ray of B{ri^i,ri) whose origin is on the bottom edge of Rec(ri^i,ri). Since 
ri-i and are two diagonal vertices of Rec{ri-i,ri), pi-i is to the right of and to the left 
of Tj-i. Further, the origin or{pi-i) is which is below r^-i and to the right of r^. 

• If rj_i is to the northeast of (see Fig. 9(3)), then if is horizontal, then the proof is 
similar to the first case; otherwise, the proof is similar to the second case. 

The lemma thus follows. □ 

Lemma 9 For any i with ^ < i < k — 1, if ri is to the southwest of ri-i, then Vi-2 is to the right 
of the rectangle Rec{ri-.i,ri) and vi is below Rec{ri-i,ri). 

Proof: Suppose is to the southwest of rj_i. We only prove that Vi-2 is to the right of the 
rectangle Rec{ri-i,ri). The case that Vi is below Rec{ri-i,ri) can be proved analogously. 

Note that Vi-2 G S(ri_2, n-i). We discuss the three possible relative positions of ri-2 and rj-i. 
By Lemma 6, rj_2 may be to the southeast, northwest, or northeast of rj_i. Since Vi is to the 
southwest of Tj-i, to prove Vi-2 is to the right of Rec{ri-i,ri), it suffices to show that Vi-2 is to 
the right of rj-i. 

• If ri-2 is to the southeast of rj_i, then by Lemma 7, Vi-2 is on the ray of i?(rj_i, rj_2) 
vertically going south, i.e., pi-2 is vertical. By Lemma 8, r^-i is to the left of Pi-2- Since 
Vi-2 G Pi-2, Vi-2 is to the right of rj_i. 

• If ri-2 is to the northwest of rj_i, then by Lemma 7, Vi-2 is to the right of Rec{ri-2,ri-i), 
and thus to the right of Rec{ri^i,ri). 

• If ri-2 is to the northeast of rj_i, then the rectangle i?ec(rj_2, Tj-i) is to the northeast of 
Rec{ri-i,ri). If Vi-2 is on BM{ri-2,ri-i), then since Vi-2 is inside Rec{ri-2,ri-i), Vi-2 is to 
the right of Rec{ri-i,ri); otherwise, the proof is similar to the above two cases. 

The lemma thus follows. □ 
Recall that when sketching the algorithm in Section 5.1, we mentioned that the origins of the 
rays involved somehow form a staircase structure. The next lemma states this important fact. 

Lemma 10 For any i with 2 < i < k — 1, or{pi-i) is to the northeast of or{pi). 
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Figure 10: Illustrating the case when rj is to the northeast of rj+i and or{pj) 7^ Vj. 

Proof: We first discuss a scenario that will be used later in this proof. Consider any two consecutive 
roots rj and r^+i in R, 1 < j < k — 1, with or{pj) ^ Vj. Then based on our discussion above, it 

must be the case that r^+i is to the southwest of rj, cd intersects the rectangle Rec{rj,rj^i), and 
or{pj) is a point on an edge of Rec{rj,rjj^i). Let Zj be the intersection of cd and the right edge of 
Rec{rj,rj^i) (see Fig. 10). The origin or{pj) can be either on the right edge or the bottom edge 
of Rec{rj,rj+i). In either case, or{pj) must be both below and to the left of Zj, i.e., zj is to the 
northeast of or{pj). 

Consider any i with 2 < i < k — 1. To prove the lemma, depending on whether or(pi_i) = ■Uj_i 
and whether or{pj) = Vi, there are four cases. 

1. If or{pi-i) = Vi-i and or{pi) = Vi, then since Vi^i and Vi are on cd in the order from c to d, 
Vi-i is to the northeast of Vi, and thus or{pi-i) is to the northeast of or{pi). 

2. If or(pi-i) = Vi-i and or{pi) ^ Vi, then by our discussion at the beginning of this proof, rj+i 
is to the southwest of r^, the rectangle i?ec(rj, rj+i) intersects cd, and the point Zi is to the 
northeast of or (pi). Further, since rj+i is to the southwest of r^, by Lemma 9, Vi-i is to the 
right of Rec{ri,ri-\.i) and thus to the right of Zi. Since Vi-i is to the right of Zi and both 
Vi-i and Zi are on cd, Vi-i is to the northeast of Zj. Therefore, or{pi-i) (= Vi-i) is to the 
northeast of or{pi). 

3. If or(/9i_i) 7^ Vi-i and or{pi) = Vi, then the analysis is somewhat similar to the second case. 

4. If or(pi_i) 7^ Vi-i and or{pi) 7^ Vi, then rj_i is to the northeast of and is to the northeast 
of Tj+i. Hence, the rectangle i?ec(r.j_i, rj) is to the northeast of i2ec(rj, rj+i). Since or(pj_i) is 
on Rec{ri-i,ri) and or{pi) is on Rec{ri, r^+i), we also obtain that or(pi_i) is to the northeast 
of or{pi). 

The lemma thus follows. □ 

Lemma 11 Consider any root ri ^ R with 1 < i < k. For any ray pj, if j < i — 1 and pj is 
vertical, then pj is to the right ofri; if j > i and pj is horizontal, then pj is below ri. 

Proof: WLOG, assume i < k. Consider the ray pi, which is on S(rj, rj_|-i). By Lemma 8, the 
origin or{pi) is below r^. By Lemma 10, for any ray pj with j > i, or(pj) is below or(pi) and thus 
is below Tj. Hence, if pj is horizontal, then pj must be below r^. 

By an analogous analysis, we can show that if j < i — 1 and pj is vertical, then pj is to the right 
of Tj. We omit the details. The lemma thus follows. □ 

Note that in any SPM, a common boundary of two adjacent cells C(r) and C(r') is a subset of 
the bisector B(r,r'). 

For any two consecutive roots r^-i and ri in R, 2 < i < k, the vertex Vi-i divides B{ri-i,ri) 
into two portions; we denote by Bi,ay{ri-i,ri) the portion that goes inside bay{cd) following Vi-i. 
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A key to building VD{bay{cd)) is to compute the interactions among all i36aj^(rj_i, rj)'s, for i = 
2,3, ... ,k, inside bay{cd). Note that if Vi-i is on a ray of B{ri-i,ri), then Bbay{ri-i,ri) is the ray 
Pi-\; otherwise, Vi-i is on BM{ri-i,ri) (i.e., the middle segment of B{ri-\,ri)), and Bhay{ri-i,ri) 
consists of a portion of BM{ri~i,ri) in Rec{ri-\,ri) (i.e., the line segment Vi-\or{pi-\)) and the ray 
Pi-i. Lemma 12 below shows that the portion of B]\{{ri^i,ri) which is inside hay{cd) will appear 
in SPM{F) (and thus in VD{bay{cd))), implying that we can simply keep it when computing 
VD{bay{cd)) and we only need to further deal with the rays for i = 1, 2, . . . , A; — 1. Thus, dealing 
with the rays pi is the main issue of our algorithm (as discussed in Section 5.1). 

Lemma 12 For any two consecutive roots rj_i and ri in R, 2 < i < k, ifvi-i lies on BM{ri-i,ri), 
then the portion of BMiri-i,ri) inside bay{cd) will appear in VD{bay{cd)). 

Proof: Consider two consecutive roots r-i^i and rj va R, 2 < i < k, with Vi-i lying on BM{T~i-i,ri). 

Denote by B'^^ the portion of BM{ri-i,ri) inside bay{cd). Recall that BM{ri-i,ri) is an open 
segment that does not contain its endpoints and is strictly inside Rec{ri-i,ri). To prove the lemma, 
it suffices to show that for any two roots rj and Vh in R with {rj,rh\ 7^ {fi-i^fi}^ if a portion of 
B{rj,rh) appears in SPM{T), then that portion does not intersect B'j^. 

By Lemma 6, may be to the southeast, or northwest, or southwest of rj_i. Since cd is 
positive-sloped, if is to the northwest or southeast of rj-i, then cd cannot intersect the rectangle 
Rec{ri-i,ri) and thus Vi-i cannot lie on BM{ri-i,ri). Therefore, the only possible case is that 
is to the southwest of Vi-i. 

First, we assume i — 1 > 2 and consider the root rj_2- We discuss the possible relative positions 
of rj_2 with respect to rj_i. Recall that the bisector portion Bi,ay{ri-2ji"i-i) either is /9j_2 or consists 
of Vi-2or{pi-2) and Pi-2- Note that in either case, when moving along Bi,ay{ri-2,ri-i) from Vi-2, 
Bbay{fi~2,ri-i) is monotonically increasing in the x-coordinates. Hence, Vi-2 is a leftmost point of 
Bbay{i'i-2-, I'i-i). Since is to the southwest of ri_i, by Lemma 9, Vi-2 is to the right of Rec{ri-i,ri) 
and thus is strictly to the right of B'j^. Hence, Bi,ay{'fi-2-,fi-i) cannot intersect B'j^j. 

For any pair of consecutive roots r^^i and rj in R, 2 < j < i — 2, similarly, when moving from 
Vj-i along Bbay{rj-i,rj), Bbay{rj-i,rj) is monotonically increasing in the x-coordinates. Since Vi-2 
is strictly to the right of B'j^ and Vj^i is to the right of Vi^2, Bbayifj-i,rj) cannot intersect B'j^. 

Let Ri = {ri,r2, . . ■ , Ti_i} and R2 = ^'i+i, • • • , r^}- (Note that since R may be a multi-set, 
Ri and R2 possibly contain the same physical root, but this is not important to our analysis.) 

For any two different pairs of consecutive roots rj-i,rj and rt-i,rt with 2 < j < i — 1 and 
2 < t < z — 1, it is possible that Bbay{rj~i,rj) and Bbay{rt-i,rt) intersect in SPM{T); if that 
happens, then let B' be the resulting bisector. It is not difficult to see that B' must be going 
in a direction between the original directions of Bi,ay{rj-i,rj) and Bbay{rt-i-,Tt)- Since neither 
Bbay{rj-i,rj) nor Bbayirt-i,rt) intersects B'j^, B' cannot intersect B'j^. We can further consider 
the possible intersection between B' and the bisector of another two roots in R\ in the manner as 
above, and show likewise that the new bisector thus resulted cannot intersect B'^. 

The above argument shows that for any two roots r j and rt in i?i such that a portion of B{r j, rf ) 
appears in VD{bay{cd)), that portion does not intersect B'f^. By a similar argument, we can also 
show that for any two roots rj and rt in R2 such that a portion of B{rj,rt) appears in VDibay{cd)), 
that portion does not intersect B'j^. 

It remains to show that for any two roots rj G Ri and rt G R2 such that {rj,rt} 7^ {ri-i,ri} 
and a portion of B{rj,rt) appears in VD{bay{cd)), that portion does not intersect B'j^^. Note that 
the case of B{rj,rt) (partially) appearing in VD{bay{cd)) can occur only after Bhayiri-i,ri) is 
"blocked" by an intersection between Bbay{ri-i,ri) and the bisector of two roots in i?i or two roots 
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Fi gure 11: Illustrating an example of BM{ri-i,ri) intersecting both cd (at Vi-i) and d (at p). The line segment 
Vi-ip divides bay{cd) into bayi and bay2. 

in i?2- Since the bisector of any two roots in Ri or any two roots in R2 cannot intersect B'j^, the 
portion of B{r.j,rt) appearing in VD{bay{cd)) cannot intersect 5^ either. 

The lemma thus follows. □ 
The observations presented above help determine the behaviors of the bisectors for the roots in 
R (e.g., the properties of the rays pi,p2, ■ ■ ■ ,Pfc-i), which are crucial to constructing VD{bay{cd)). 
They form a basis for both showing the correctness and the efficiency of our algorithm in Section 
5.3. For example, Lemma 10 can help conduct a set of ray shooting operations in linear time, and 
Lemma 12 allows us to decompose the problem into certain subproblems with good properties. 

5.3 The Algorithm for Computing VD{hay{cd)) 

In this subsection, we present our algorithm for computing VD{bay{cd)), i.e., computing the Voronoi 
region VD{r) for each root r £ R. 

As shown in [29, 30], a key property of the problem in the Li metric is: There exists an SPM 
such that each edge of the SPM is horizontal, or vertical, or of a slope 1 or —1. As shown below, 
the curves involved in specifying VD{hay{cd)) consist of only line segments of slopes 0, 00, and — 1 
(there is no +1, which is due to the assumption that cd is positive-sloped). A line (segment) is 
said to be (—l)-sloped if its slope is —1. Our algorithm needs to perform some vertical, horizontal, 
and (— l)-sloped ray shooting queries, whose total number is 0{k). By exploiting some properties 
of our problem shown in Section 5.2, we conduct all ray shootings in a global manner in totally 
0{n' + k) time. 

The pseudo-code of Algorithm 1 summarizes the entire algorithm. 

Before describing the main algorithm, we discuss some preprocessing work as well as some basic 
algorithmic methods that will be used later in the main algorithm. 

5.3.1 Preliminaries and Preprocessing 

By Lemma 12, for any two consecutive roots rj_i and ri in R, 2 < i < k, if the middle segment 
BM{ri-i,ri) of their bisector intersects cd (at Vj-i), then we can "separately" process the portion 
of BM{ri-i,ri) inside Rec{ri-i,ri), as follows. Let d be the boundary of bay{cd) minus cd, i.e., d 
consists of all edges of bay{cd) except cd. 

Clearly, Vi-i divides BMifi-i^ri) into two portions; one portion does not contain any point 
in bay{cd) and the other contains some points in bay{cd). Denote by B'j^(ri-i,ri) the portion 
that contains some points in bay{cd). Thus, B'j^{ri-i,ri) is a line segment and Vi-i is one of its 
endpoints (and or{pi-i) is the other endpoint). We first determine whether B'j^{ri-i,ri) intersects 
d, by performing a — 1-sloped ray shooting operation. Specifically, we shoot a ray p originating 
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Algorithm 1: Computing a shortest path map for bay{cd) 



Input: hay{cd), R = {ri, r2, . . . , r^}, and SPM{A4) vertices vi,V2, ■ ■ ■ , v^-i- 
Output: A shortest path map on bay(cd) with respect to the source point s. 



/* Preprocessing 

Compute the ray set ^ = {pi , p2 , • 



*/ 



for each 1 < i < 



2 Compute the Une segment ViOr{pi) for each 1 < i < k — 1 if Vi ^ or{pi) ; 

3 Compute the horizontal visibihty map HM{bay{cd)) and the vertical visibility map VM{bay{cd)) ; 

4 Compute the trapezoid in HM{bay{cd)) that contains or{pi) for each 1 < i < A; — 1 ; 
/* The main algorithm 

5 p* ■<— c, 5 0, Q {pi,p2, • ■ ■ ,pk-i} ; /* Q is a queue storing the rays. 

6 while Q is not empty do 
Consider the first ray p in Q and remove it from Q ; /* Assume p is on B[rj,ri) with i > j . 
if p is vertical then 

I Push p onto the top of S, and exit the current loop ; 
else /* p ±s horizontal. 

Compute the target point tp{p) ; 
if S is empty then 

The Voronoi region VD{rj) is determined with or{p)tp{p) ; 
p* -ir- tp{p), and exit the current loop ; 
else /* S is not empty; assume p' C B{rt,rj) with j > t is the ray at the top of 



S. 



*/ 
*/ 

*/ 



*/ 



*/ 



*/ 



assume 

Scan d{p* ,tp{p)) to compute the target points on d{p* ,tp{p)) of the rays in S ; 
if tp{p') is before tp{p) (i.e., tp{p') has been computed) then 

Determine the Voronoi regions for the roots defining the rays in S ; 
Pop all rays out of S ; 
p* <— tp{p), and exit the current loop ; 
else/* tp{p') is not before tp{p) (i.e., tp{p') has not been computed). 
Determine the Voronoi region VD{rj) ; 

Let p be the intersection of p and p' , and q be the intersection of the horizontal line through 
ri and the vertical line through r*; let p' be the other intersection of BM{rt,ri) and the 
boundary of Rec{p, q) than p; 

Move from p along pp' in HM{bay{cd)) until either p' or d is encountered first; 
if d is encountered ( say, at the point z ) then 

Scan d{tp{p),z) to compute the target points on d{tp{p),z) of the rays in S ; 
Determine the Voronoi regions for the roots defining the rays in S ; 
Pop all rays out of S ; 
p* <— z, and exit the current loop ; 
else /* d is not encountered. */ 
Pop p' out of S ; 

if p' is on the bottom edge of Rec{p, q) then 

I Push the ray originating at p' and going south onto the top of S ; 
else /* p' is on the right edge of Rec{p, q) . */ 
|_ Add the ray originating at p and going east to the front of Q ; 

p* *r- tp{p), and exit the current loop ; 



37 For each n £ compute the SPM on the Voronoi region VD{ri) with respect to n 
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Figure 12: Illustrating the horizontal visibility map of a simple polygon. 

at Vi-i and passing through the other endpoint of B'j^(ri-\,ri). If the length of the portion of p 
between Vi-\ and the first point p on d hit by p is larger than the length of B'j^{ri-i,ri), then 
B'j^{ri-\,ri) does not intersect 5, and we do nothing. Otherwise, B'j^{ri-i,ri) intersects d (at 
the point p). By Lemma 12, the line segment Vi^ip appears in SPM{T). Also, Vi-\p partitions 
bay{cd) into two simple polygons (see Fig. 11); one polygon contains cvi^i as an edge, which 
we denote as bayi, and we denote the other polygon as 5aj/2- Let Ri = {ri, r2, . . . , rj-i} and 
Ri = {'r'i,ri-\-i, . . . ,rk}. (Note that since R may be a multi-set, Ri and R2 possibly refers to the 
same physical root, but this is not important to our algorithm.) Since Vi-ip is in VD(bay(cd)), 
it is not difficult to see that for any point q in bayi, there is a root r G Ri such that a shortest 
path from s to g goes through r. Similarly, for any point q in bay2, there is a root r £ R2 such 
that a shortest s-q path goes through r. This implies that we can divide the original problem of 
computing VD{bay{cd)) on bay{cd) and R into two subproblems of computing VD(6ayi) on bayi 
and Ri and computing VD (601/2) on 607/2 and i?2- 

If we process each pair of consecutive roots in R as above, then the original problem may 
be divided into multiple subproblems, each of which has the following property: For any pair of 
consecutive roots rj_i and in the corresponding root subset of R, if BM{ri-i,ri) intersects cd, 
then B'j^{ri-i,r,i) does not intersect d and is contained in the corresponding subpolygon of bay(cd); 
further, B'j^j{ri^i,ri) is in VD{bay{cd)) and has been computed. 

To perform the above process, a key is to derive an efficient method for the — 1-sloped ray 
shooting operations. For this, we choose to check all pairs of consecutive roots in R in the order 
of ri,r2, . . . ,rjk. In this way, it is easy to see that the ray shootings are conducted such that the 
origins of the rays are sorted along cd from c to d. This is summarized by the next observation. 

Observation 5 The preprocessing conducts 0{k) —1-sloped ray shooting operations that are orga- 
nized such that the origins of all rays are on cd ordered from c to d. 

We show next that the ray shootings for Observation 5 can be done in 0{n' -\-k) time. Since the 
origins of all rays in Observation 5 are sorted on cd, we can perform the ray shootings by computing 
the visible region of bay{cd) from cd along the direction of these rays. This can be easily done by 
a visibility algorithm on a simple polygon (e.g., [1, 23, 28]). Below, we give a different algorithm 
for a more general problem; this more general result is needed by the main algorithm. 

Given a simple polygon P, the horizontal visibility map of P contains a horizontal line segment 
inside P through each vertex of P, extending as long as possible without properly crossing the 
boundary of P (such line segments are called the diagonals; see Fig. 12). The vertical visibility map 
with vertical diagonals is defined similarly. Each region in a visibility map is a trapezoid (a triangle 
is a special trapezoid). A visibility map of a simple polygon can be computed in linear time [3]. 

For a ray p with its origin in bay{cd) (inside it or on the boundary), the boundary point of 
bay{cd) that is not the origin or{p) hit by p first is called the target point of p, denoted by tp{p). 
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Recall that d is the boundary of hay{cd) excluding the edge cd. In the rest of this paper, unless 
otherwise stated, a ray in our discussion always has its origin in hay{cd) and its target point on d. 

We say that m parallel rays p'^^ p'^-, ■ ■ ■ ■, p'm ^'^^ target-sorted if we move from c to d (clockwise) 
on 5, we encounter the target points of these rays on d in the order of tp{p[), tp{p2), . . . , tp{p'^). 

Given a set of m target-sorted parallel rays p'l, P2, ■ ■ ■ , p'm for bay{cd) whose origins are in bay{cd) 
and whose target points arc on d, below we present a visibility map based approach for computing 
their target points in 0{n' + m) time (recall that n' is the number of vertices of bay{cd)). 

WLOG, we assume that the rays are all horizontal. We first compute the horizontal visibility 
map of bay{cd) in 0{n') time. Then, starting from the vertex c, wc scan d and check each edge 
e of 9 and the trapezoid t(e) of the visibility map bounded by e, to sec whether the next ray p[ 
(initially i = 1) is in the trapezoid t{e) and can hit the edge e. Once the target point of the ray 
is found, we continue with the next ray p'^^i- Clearly, the time for computing all target points is 
0{n' + m). Thus, we have the following result. 

Lemma 13 Given a set of m target-sorted parallel rays for bay{cd) whose origins are in bay{cd) 
and whose target points are on d, their target points can be computed in 0{n' + m) time. 

For the ray shootings in Observation 5, it is easy to see that these rays are target-sorted. Thus, 
by Lemma 13, their target points can be computed in 0{n' + k) time (of course, these ray shootings 
can be done by using the visibility algorithms in [1, 23, 28], which do not compute a visibility map). 
We present the above visibility map based technique because our main algorithm in Section 5.3.2 
will need it. 

In addition, as part of the preprocessing for our main algorithm, we also compute the horizontal 
visibility map HM{bay{cd)) and the vertical visibility map VM{bay{cd)) of bay{cd). Further, for 
each 1 < i < A; — 1, we compute the trapezoid of the horizontal visibility map HM{hay{cd)) that 
contains the origin or{pi) of the ray pi, in totally 0{n' + k) time, in the following way. 

Recall that or{pi) is either Vi or in the interior of bay{cd). In the latter case, or{pi) is an 
endpoint of the line segment i3^(ri,ri+i) = Vior{pi) whose slope is —1, and the position of or(pi) 
has been determined earlier by the — 1-sloped ray shooting operations. By Lemma 10, all origins 
or{pi),or{p2), . . . , or(/9fe-i) are ordered from northeast to southwest. Further, or^piYs are all visible 
from cd along the direction of slope —1. Thus, it is not difficult to show that if we visit the trapezoids 
of HM(bay{cd)) by scanning the edges of d from c to d and looking at the trapezoids bounded 
by each edge, then the trapezoids containing such or(pi)'s are encountered in the same order as 
or{pi),or{p2), ■ ■ ■ ,or{pk-i)- This implies that we can use a similar algorithm as for computing 
the target points of target-sorted parallel rays on d (i.e., scanning d from c to d and checking the 
trapezoids of HM(bay{cd)) thus visited along d) to find all the sought trapezoids, in 0{n' -\-k) time. 

The above discussion leads to the following lemma. 

Lemma 14 The preprocessing on bay{cd) takes 0{n' -\- k) time. 

In the main algorithm, the horizontal visibility map HM{bay{cd)) will be used to guide the 
main process. More specifically, during the algorithm, we traverse inside bay{cd) following certain 
rays, and use HM{bay{cd)) to keep track of where we are (i.e., which trapezoid of HM{bay{cd)) 
contains our current position). The vertical visibility map VM{bay{cd)) will be used to compute 
the target points of some target-sorted vertical rays using the above visibility map based approach. 

For any two points a and b on d with a lying on the portion of d from c clockwise to b, we 
denote by d{a, b) the portion of d between a and b and say that a is before 6 or 6 is after a. 
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5.3.2 The Main Algorithm 

After the preprocessing, the problem of computing VD(bay{cd)) with the root set R may be di- 
vided into multiple subproblems and we need to solve each subproblem. For convenience of the 
forthcoming discussion, we assume that the original problem on bay{cd) with R is merely one 
such subproblem, i.e., for any two consecutive roots and rj+i in R, if Vi G BM{r^,ri-^-l), then 
B']^{ri,ri+i) (= Vior{pi)) lies completely in VD{bay{cd)) and has been computed. Recall that 
in the preprocessing, we have already computed the trapezoid of the horizontal visibility map 
HM{bay{cd)) that contains the origin or{pi) of the ray pi, for each 1 < i < k — 1. Observation 6 
below summarizes these facts. 

Observation 6 After the preprocessing, 

• for any two consecutive roots Vi andri-\-\ in R, ifvi G -BM(ri, rj+i), then their bisector portion 
B^(rj,rj+i) (= ViOr{pi)) has been computed; 

• for each 1 < i < k — 1, the trapezoid of HM(bay{cd)) that contains the origin or{pi) of the 
ray pi is known. 

As discussed before, our task is to handle the interactions among the rays pi for all i = 
l,2,...,fc-l. 

In the algorithm, we need to compute the target points for 0{k) horizontal and vertical rays. 
The main procedure is guided by the horizontal visibility map HM{bay{cd)) so that the target 
point of each horizontal ray can be determined in constant time. For the vertical ray shootings, 
we use the visibility map based approach with the vertical visibility map VM{bay(cd)). Note that 
the vertical ray shootings will occur in an online fashion in the algorithm. We will show that the 
vertical rays involved are target-sorted. To compute the target points for these vertical rays, the 
algorithm maintains a reference point, denoted hy p* . Initially, p* = c. Then during the algorithm, 
p* will be moved forward along d from c to d, i.e., every time p* is moved on d, its new position is 
always after its previous position. In this way, the target points of all vertical rays are computed 
in totally 0(n' + k) time (recall that n' is the number of obstacle vertices of bay{cd)). 

Let ^' = {/Oi, /92) • • • ) Pk-i}- We process the rays of ^' incrementally in the order of pi, P2, ■ ■ ■ , Pk-ij 
whose origins are ordered from northeast to southwest by Lemma 5 10. By Observation 4, each ray 
in is either horizontally going east or vertically going south. We say that initially all rays are 
active and the entire bay{cd) is active. In general, the active rays are used to decompose the active 
region of bay{cd). During the algorithm, some portion of bay[cd) will be implicitly set as inactive, 
which means that each point of such a region is in the Voronoi region of a root that has been de- 
termined. The active region of bay{cd) at any moment of the algorithm always forms a connected 
simple polygon, a fact that we will not explicitly argue in the following algorithm description. Sim- 
ilarly, some rays will be set as inactive, meaning that they will no longer be involved in the further 
decomposition of the current active region of bay{cd). When the algorithm terminates, the entire 
bayicd) is inactive and all rays of ^ are inactive. Note that setting a region or a ray as inactive is 
done implicitly and is used only for our analysis. Since each ray in ^ lies on the bisector of two 
roots in R, we say that the two roots define the ray. 

We start with the first ray pi. If pi is horizontal (going east), then since or(pi) is the most 
northeast origin, no other ray in can intersect it. Let p be the target point of pi on d (see 
Fig. 13). Clearly, p can be found in 0(1) time since we already know the trapezoid in HM{bay{cd)) 
that contains or{pi) by Observation 6. Denote by a the portion of B(ri,r-2) between vi and p. 



28 



Figure 13: Illustrating an example of pi being horizontal. 

Note that a is either the hne segment vip (if vi = or(pi)), or the concatenation of the two hne 
segments vior{pi) and or{p\)p. In either case, a partitions bay{cd) into two simple polygons. One 
of them contains cuT as an edge and we denote it by bayi (see Fig. 13). We claim that bayi is 
the Voronoi region of ri, i.e., VD{ri) = bayi. Indeed, by the above analysis and Lemma 12, a is 
in VD{bay{cd)), implying that for any point q G bayi, there is a shortest path from s to q via ri. 
The claim thus follows, and VD{ri) is determined. Wc then set the ray pi and the region of bayi 
as inactive. Hence, the active region of bay{cd) becomes bay{cd) \ bayi, which needs to be further 
decomposed. In addition, we move the reference point p* from c to p (= tp{pi)). We then continue 
with the next ray p2- 

If pi is vertical (going south), then we push pi onto a stack 5* (initially, iS = 0), and let the 
reference point p* stay at c. We then continue with the next ray p2- 

We will show below that our algorithm maintains the following general invariants, which are 
used to prove the correctness of the algorithm. Suppose the current moment of the algorithm is 
right before the next ray p is considered, and assume p lying on the bisector B{rj,ri) with i > j. 
The stack S may be non-empty; if 5" = 0, then the invariants below related to any rays in S are 
not applicable. Let p' be the ray at the top of S, and suppose p' lies on B{rt,rti) with t' > t. 

Invariant Properties: (1) All rays in S are active and vertically going south. (2) The origins 
of all rays in S from top to bottom are ordered from southwest to northeast. (3) The origin of the 
next ray to be considered by the algorithm (i.e., p) is to the southwest of the origin of the ray at 
the top of S (i.e., p'). (4) The two indices j = t' . (5) For each ray p" in 5 U {p}, suppose p" lies on 
the bisector B{rj',rii) of two roots rji and r^i with i' > f; then the portion of the boundary of the 
Voronoi region yD(rj') (resp., VD{rji)) from Vi'^i (resp., vji) to the origin or{p") of p" has already 
been computed. (6) For each ray p" in S, suppose it lies on the bisector B{rji,rii) of two roots r^/ 
and Tj' with i' > f; then r^/ is to the right of p" and rj/ is to the left of p" . (7) The root rt is to the 
left of all rays in 5 \ {p'} (recall that p' C B{rt,rti) with t' > t). (8) For any two consecutive rays 
p[ and P2 in S such that p'l is closer to the top of S, suppose p'^ is on B{ri^,ri^) for i2 > ii and 
p'2 is on B{rj-^, rj,^) for j2 > ji; then ii = j2- (9) The target points of all rays in 5* from bottom to 
top are ordered clockwise on d (i.e., from c to d). (10) If p is vertical, then the target point tp{p) 
of p is after the target point tp{p') of p' on d. (11) If the target point of any ray in S has not been 
computed yet, then the target point of that ray is after the reference point p* (i.e., on d{p*,d)). 
(12) The target point tp{p) is after p*. (13) Suppose p" is the first horizontal ray in that will 
be considered by the algorithm in a future time from now; then its target point tp{p") is after p* . 
(14) The trapezoid in the horizontal visibility map HM{bay{cd)) that contains the origin or{p) of 
the ray p is known. 

Now consider the moment that is right after we finish processing pi and before we consider 
P2- Based on the processing of pi discussed above, either pi is horizontal and S is empty, or pi is 
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vertical and S = {pi}. Lemma 15 below shows that in either case, all invariants of the algorithm 
hold. We intend to use the proof of Lemma 15 as a "warm-up" for the analysis of the more general 
situations later. Observation 7 follows from the definitions of the rays in * and Lemma 10. 

Observation 7 The target points of all rays in ^ are on d. For any two rays rj and ri in ^ with 
i > j, if rj is horizontal or ri is vertical, then tp{ri) is after tp{rj) on d. 

Lemma 15 At the moment after pi has been processed and before p2 is considered, all invariants 
of the algorithm hold. 

Proof: Recall that pi is on B{ri,r2) and p2 is on B{r2,rz), and the reference point p* is at the 
target point tp{pi) if pi is horizontal and at the vertex c otherwise. 

We first discuss the case when pi is horizontal, in which S is empty and p* = tp(pi). Invariants 
(1) through (11) except (5) simply follow since they are all related to some rays in S. For Invariant 
(5), we only need to consider p2 C i?(r2, ra), i.e., we need to show that the portion of the boundary 
of the Voronoi region VD{r3) from V2 to the origin or{p2) of p2, which is also the boundary portion 
of the Voronoi region VD{r2) from V2 to or{p2), has already been computed. Denote this boundary 
portion by a. Note that a is the portion of B(r2,r3) between V2 and or(p2)- Recall that or{p2) is 
either V2 or not. If or{p2) = V2, then we are done since a is just a single point V2- Otherwise, V2 must 
be on -Bm(?'2,?'3) and a is B'j^{r2,r3) (= V2or{p2)), which has been computed in our preprocessing 
by Observation 6. Hence, Invariant (5) follows. 

For Invariant (12), we need to show that tp{p2) is after p* = tp{pi), which is true due to 
Observation 7 and pi being horizontal. For Invariant (13), let i > 1 be the smallest index such that 
Pi G * is horizontal. If there is no such i, then Invariant (13) trivially holds; otherwise, we need to 
prove that tp{pi) is after p* = tp{pi), which is true due to Observation 7 and pi being horizontal. 
For Invariant (14), we need to show that the trapezoid of HM{bay{cd)) containing or(p2) is known, 
which is true by Observation 6. Hence, when p\ is horizontal, all invariants hold. 

We then discuss the case when pi is vertical, in which S = {pi} and the reference point p* = c. 
Invariants (1) and (2) simply follow. By Lemma 10, or{p2) is to the southwest of or(pi), and thus 
Invariant (3) holds. Invariant (4) is obvious. For Invariant (5), we need to consider both pi and p2. 
The proof is similar to that for the case when pi is horizontal, and we omit it. For Invariant (6), we 
need to show that ri is to the right of pi and r2 is to the left of pi, which is true due to Lemma 8 
and pi being vertical. Invariant (7) simply follows since pi is the only ray in S. Invariants (8) and 
(9) trivially hold since S has only one ray. For Invariant (10), we need to show that if p2 is vertical, 
then tp{p2) is after tp(pi), which is true by Observation 7. For Invariant (11), note that the target 
point tp{p\) has not been computed. Since p* = c. Invariant (11) trivially holds. Invariants (12) 
and (13) also easily hold since p* = c and the target points of all rays in * are on d. For Invariant 
(14), we need to show that the trapezoid of HM{bay{cd)) that contains or{p2) is known, which is 
true by Observation 6. 

We hence conclude that all invariants of the algorithm hold. □ 
As an implementation detail, although we view 5 as a stack, we represent 5 as a doubly-linked 

list so that we can access the rays in S from both the top and the bottom of S. But, we always 

pop a ray out of S from its top and push a ray onto S at its top. Next, we discuss the general 

situations of our algorithm. 

Suppose our algorithm just starts to process a ray pj G 4', i > 1, which lies on the bisector 

B{ri, Tj+i), and all invariants of the algorithm hold right before pi is processed. There are a number 
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of cases and subcases to consider, depending on whether pi is vertical or horizontal, whether S is 
empty, and the intersecting consequences between pi and the rays in S (if S" 7^ 0), etc. 

Case 1: is vertical (going south). Then we simply push pi onto the top of S and the reference 
point p* is not changed. The algorithm then continues with the next ray pi+i G * in this situation. 
Lemma 16 below shows that all invariants of the algorithm hold. 

Lemma 16 // the ray pi is vertical, then at the moment after pi is processed and before pi+i 
is considered, all invariants of the algorithm hold. 

Proof: Note that pi+i is on i?(ri_|_i, ri+2). Let ^ be the moment right after pi is processed and ^' 
be the moment right before pi is considered. Thus, from ^' to ^, the only change to S is that we 
push Pi onto the top of S. The proof below is based on the assumption that S has at least two rays 
at the moment ^ (i.e., pi and at least another ray), since otherwise the invariants related to other 
rays in S than pi trivially hold. This also implies that S is not empty at the moment Let p be 
the ray at the top of S at the moment and assume p lying on the bisector B{rj,rji) with / > j. 
Invariant (1) holds since pi is vertical. 

For Invariant (2), since all invariants of the algorithm hold at the moment it suffices to show 
that or{pi) is to the southwest of or{p). Note that at the moment pi is the next ray to be 
considered by the algorithm. Thus, by Invariant (3) at the moment or{pi) is to the southwest 
of or{p). Invariant (2) thus follows. 

For Invariant (3), Lemma 10 implies that or{pi+i) is to the southwest of or{pi). 

Invariant (4) trivially holds since pi C i?(rj,rj+i) and pj+i C -B(rj+i, rj-|-2)- 

For Invariant (5), it suffices to consider the ray pi+i, i.e., to show that the portion of the 
boundary of VD(ri+i) between Vi+i and or(pi+i), which is also the portion of the boundary of 
VD{ri^2) between Vj+i and or{pi+i), has already been computed. (Note that the case for the ray 
Pi trivially holds due to Invariant (5) at the moment ^' when pi is the next ray to be considered.) 
Recall that the above boundary portion is a single point Wj+i if Vi^i = or(/7j_|_i) and is the line 
segment fi+ior(/9i+i) otherwise. By Observation 6, if fj+i 7^ or(/9i+i), then f j_|_ior(pi+i) has 
already been computed in the preprocessing. Thus, Invariant (5) follows. 

For Invariant (6), it suffices to prove that is to the right of pi and rj+i is to the left of pi, 
which follows from Lemma 8 since pi is vertical. 

For Invariant (7), it suffices to show that is to the left of p. Recall that p is on B{rj,rjf) with 
/ > j. At the moment by Invariant (4), f = i; by Invariant (6), (= rji) is to the left of p. 
Thus, Invariant (7) follows. 

For Invariant (8), it suffices to show j' = i, which has been proved above for Invariant (7). 

For Invariant (9), it suffices to show that tp{pi) is after tp{p) on d. Since pi is vertical, at the 
moment by Invariant (10), tp(pi) is after tp{p) on d. Invariant (9) thus follows. 

For Invariant (10), we need to prove that if pi+i is vertical, then tp(pj_|_i) is after tp{pi) on d, 
which follows from Observation 7. 

For Invariant (11), note that the target point tp{pi) has not been computed. We need to show 
that tp{pi) is after p* on d. At the moment by Invariant (12), tp{pi) is after p*. Further, p* has 
not been moved since the moment Invariant (11) thus follows. 

For Invariant (12), we need to show that tp{pi+i) is after p*. If pj+i is vertical, then by 
Observation 7, tp(pi^i) is after tp{pi) on d, and we have also shown above that tp{pi) is after p*; 
thus tp{pi^i) is after p*. If pj+i is horizontal, then at the moment since pi is vertical, the first 
horizontal ray in 'if to be considered by the algorithm in future is Pi+i] thus by Invariant (13), 
tp{pi+i) is after p*. Invariant (12) then follows. 
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Invariant (13) trivially holds since pi is vertical. More specifically, suppose the first horizontal 
ray in ^! that will be considered by the algorithm after the moment ^' is pj. Note that j > i. Since 
all invariants of the algorithm hold at the moment by Invariant (13), tp{rj) is after p* on d. 
Then at the moment ^, since pi is vertical, the first horizontal ray in ^ to be considered by the 
algorithm is still pj. Proving that Invariant (13) holds at the moment ^ is to prove that tp{rj) is 
after p* , which has been proved above since p* has not been moved since the moment 

For Invariant (14), we need to show that the trapezoid of HM(bay(cd)) that contains or(pj_|_i) 
is known, which is true by Observation 6. 

We conclude that all invariants of the algorithm hold at the moment ^. □ 

Case 2: pi is horizontal (going east). Let p = tp{pi). We claim that we can find p in constant time. 
Indeed, since pi is the next ray considered by the algorithm, by Invariant (14), the trapezoid of 
HM{bay{cd)) that contains or{pi) is known. The claim then follows since p is on the boundary of 
the above trapezoid. Since pi is horizontal, by Invariant (13) (at the moment right before processing 
Pi), p is after the reference point p*. Depending on whether the stack S is empty, there are two 
subcases to consider. 

Subcase 2(a): S" = 0. Then no ray in S intersects pi before it hits d (and thus no ray shooting 
for any ray pj G ^ with j < i intersects pi before hitting d). Also, for each ray pj G ^' with j > i, 
since or{pj) is to the southwest of or{pi) and pi is horizontal, pj cannot intersect pi. Hence, the 
portion or{pi)p of the ray pi appears in VD{hay{cd)). Recall that pi is on S(ri,ri+i). The portion 
of -B(ri, rj+i) between Vi and p divides the current active region of bay{cd) into two simple polygons; 
one of them contains Vi-iVi and we denote it by bayi. Further, each point in bayi has a shortest 
path to s via r^. Thus, bayi is the Voronoi region VD{ri). We then set pi and the region bayi as 
inactive. In addition, we move p* to p. We then consider the next ray p^+i . We prove below that 
all invariants of the algorithm hold right after processing pj. 

Since S is empty. Invariants (1) to (11) except (5) simply hold since they are all related to some 
rays in S. For Invariant (5), we only need to consider pi+i, which also holds by Observation 6 (the 
analysis is similar as before). For Invariant (12), we need to show that ip(pi+i) is after p* . Since 
Pi is horizontal, by Observation 7, tp{pi+i) is after tp{pi) (= p*). Thus, Invariant (12) follows. 
For Invariant (13), suppose pj G ^' is the first horizontal ray to be considered by the algorithm. 
Note that it must be j > i. We need to show that tp{pj) is after p* (= tp{pi)), which is true 
by Observation 7 since pi is horizontal. For Invariant (14), we need to show that the trapezoid 
of HM{bay{cd)) that contains or(pi+i) is known, which is true by Observation 6. Therefore, all 
invariants of the algorithm hold right after processing pi. 

Subcase 2(b): S ^ $. Then for the rays in S whose target points lie on d{p*,p), we compute 
their target points by scanning d{p*,p) from p* to p (= tp{pi)); this scanning process uses the 
visibility map based approach with VM(bay{cd)), as described in the preprocessing. By Invariant 
(9), such vertical rays (from bottom to top in S) are target-sorted. Thus, the scanning procedure 
takes linear time in terms of the number of edges of d{p*,p) and the number of target points found 
in this process. Note that the scanning procedure stops when we encounter the point p. This also 
implies that the target points of some rays in S (e.g., the ray at the top of S) are not yet found if 
they are on d after p. 

Let p be the ray at the top of S (e.g., if pi-i is vertical, then p is pi-i). Suppose p is on the 

bisector B{rj,rj') with f > j. Then right before pi is processed, by Invariant (4), i = f since 
Pi C B{ri,ri+i); by Invariant (3), or{pi) is to the southwest of or{p). Depending on whether the 
target point tp{p) is before p, there are two subcases. 
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Fi gure 14: Illustrating an example that the target points of all rays in S are before p = tp{pi). All vertical rays are 
in S. The ray p is at the top of S and p' is at the bottom of S. 

Subcase 2(b.l): The target point tp{p) is before p = tp{pi). In this case, the scanning procedure 
has found tp{p) on d{p* ,p). Then by Invariants (9) and (11), the target points of all rays in 5" have 
been obtained and all such target points are before p on d. Since p is the target point of pi, the 
above implies that all rays in S hit d before they intersect pi (see Fig. 14). Further, since all other 
active rays in ^, i.e., pi+i,pi-\-2, ■ ■ ■ ,Pk-i, have their origins to the southwest of or{pi), no ray in S 
can intersect these active rays before it hits d. This means that for each ray p' in S, the portion of 
p' between or{p') and tp{p') appears in VD(bay{cd)). Based on the discussion above, we perform a 
splitting procedure on S, as follows. 

Let p' be the ray at the bottom of S and z = tp{p') (see Fig. 14). Suppose p' is on the bisector 
B(rt,rt') with t' > t. By Invariant (5), the boundary portion of VD{rt) between vt and or{p') has 
been computed. The concatenation of the segment or{p')z and this boundary portion of VD{rt) 
splits the current active region of bay{cd) into two simple polygons. One of them contains vt-ivt 
as an edge; further, each point in this polygon has a shortest path to s via rt- Thus, the polygon 
containing vt-\vt is the Voronoi region VD{rt)- We also set the region VD{rt) as inactive. 

We then continue to process the second bottom ray in 5, in the similar fashion. This splitting 
procedure stops once all rays in S are processed. In addition, we set all rays in S as inactive and 
pop them out of S {S then becomes empty). Finally, we move the reference point p* to p (= tp{pi)), 
and consider the next ray Pi+i- By the same analysis as that for the subcase 2(a) when S is empty, 
we can prove that all invariants of the algorithm hold. We omit the details. 

Subcase 2(b.2): The target point tp{p) is not before p = tp{pi). In this case, tp{p) has not been 
found on d{p*,p) by the scanning procedure. Then, it is easy to see that p intersects pi before it 
hits d (and so may some other rays in S) . We need to consider the consequences of the intersections 
of such rays in S with pi. Recall that pi is on S(rj,ri+i) and p is on B{rj,ri) with i > j. Below 
we show how to determine the Voronoi region VD(rj) and the portion of the bisector S(rj,rj+i) in 
VD(bay{cd)). Let pi be the intersection point of pi and p (see Fig. 15). 

First of all, we determine the Voronoi region VD{ri) (see Fig. 15). Since p is the leftmost ray in 
S by Invariant (2), both the line segments or{pi)pi and or{p)pi appear in VD{bay{cd)). Since the 
ray p is in the stack S, by Invariant (5), the boundary portion of VD{ri) between Vi-i and or(p) 
has been computed, which we denote by a. At the moment right before pi is processed, since pi 
is the next ray to be considered, also by Invariant (5), the boundary portion of VD{ri) between Vi 
and or{pi) has been computed, which we denote by (3 (i.e., (i = Vi \i vi = or{pi) and j3 = Vior{pi) 
otherwise). As argued similarly in the earlier analysis, VD{ri) is the region bounded clockwise by 
a, the segment or{p)pi, the segment or{pi)pi, P, and Vi-iVi. This region of VD(rj) is then set as 
inactive. 
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Figure 15: Illustrating an example that the ray p at the top of S intersects pi (at pi) before p hits d. 

Second, we determine the portion of the bisector B(rj,ri^i) that appears in VD(cd). Since pi 
is horizontal, by Lemma 8, the root rj+i is below pi. Since p C B{rj,ri) with i > j, by Invariant 
(6), Tj is to the right of p. Therefore, the intersection point (denoted by qi) of the horizontal hne 
through rj+i and the vertical line through rj is to the southeast of pi (see Fig. 15). We first discuss 
the portion of B{rj,ri^i) contained in the rectangle Rec{pi, qi). 

Obviously, Rec{pi, qi) is contained in the rectangle i?ec(rj, rj+i). Thus, the portion oi B{rj,ri^i) 
in Rec{pi,qi) is a portion of the middle segment of B{rj,ri^i). Further, since pi is the intersection 
of Pi and p, pi is at the intersection of B{ri,riJ^i) and B{rj,ri). Thus, pi is on B{rj,ri^i). 

We claim that rj+i is to the southwest of rj. This can be proved by showing that rj+i is to 
the southwest of pi and pi is to the southwest of rj. Indeed, since rj+i is below pi and cd is 
positive-sloped, pi must be to the right of rj+i, which also implies that rj+i is to the southwest of 
pi . Similarly, we can show that pi is to the southwest of rj . 

Because rj+i is to the southwest of rj, the middle segment of B{rj,ri^i) is — 1-sloped. Denote 
by B']y.j{rj,ri^i) the portion of i?(rj,rj+i) contained in Rec{pi,qi). Based on the above analysis, 
B'j^,j{rj,ri^i) is a —1-sloped line segment with an endpoint at pi and the other endpoint on one of 
the two edges of Rec{pi, qi) incident to qi (see Fig. 15). Below we prove that B';^_j{rj , ri+i) ribay{cd) 
(i.e., the portion of B'j^j{rj,ri+i) contained in bay{cd)) appears in VD{bay{cd)), implying that we 
should keep this portion of i?^^(rj, rj+i). The proof is similar to that for Lemma 12 and hence 
we only sketch it here. For convenience, we view B'j^j{rj,ri^i) as the open segment that does not 
contain its two endpoints. 

It suffices to show that B'j^j{rj,ri^i) does not intersect any current active ray. Consider any 
current active ray p' , p' {p, pi}. Then p' either is in or is a ray pt ^ with t > i. 

• If € S", then p' is vertical by Invariant (1). By Invariant (7), p' is to the right of the root r^, 
and thus to the right of the rectangle Rec{pi,qi). Hence, pt does not intersect i3^(rj,rj+i) 
since S^^(rj, rj+i) is strictly inside Rec{pi,qi). 

• p' = Pt ^ with t > i, then there are two subcases. 

If Pt is horizontal, then by Lemma 11, pt is below rj+i, and is thus below the rectangle 
Rec{pi,qi). Hence, pt does not intersect i?^j(rj, rj+i). 

If Pt is vertical, then by Lemma 10, the origin or(pi) is to the northeast of or{pt). Clearly, 
or{pi) is to the left of Rec{jpi,qi) and thus or{pt) is to the left of Rec{pi,qi). Since pt is 
vertical, pt is also to the left of Rec{pi^qi). Hence, pt does not intersect S^^(rj, rj+i). 
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Figure 16: Illustrating an example that _BJ^f (rj, r^+i) {= pip'i) intersects d (first at z). 

The above argument shows that all active rays cannot intersect B'pj{rj,ri^i). Hence, the portion 
of S^j(rj,rj+i) contained in bay{cd) must appear in VD{hay{cd)). 

Then, we compute B'j^j{rj,riJ^i) in 0(1) time, and let Bpj{rj,ri+i) = pip'i (see Fig. 15). Note 
that p[ is either on the right edge or the bottom edge of Rec{pi,qi). 

However, pip'i may intersect d. To determine whether such intersection occurs, we move among 
the trapezoids in the horizontal visibility map HM{hay{cd)) from the endpoint pi of B'j^j{rj,ri^i) 
along the segment pip'i, as follows. 

Note that the portion of the ray pi between its origin or(pi) and its target point p = tp{pi) 
is contained in a single trapezoid of HM{bay{cd)), i.e., the trapezoid containing or{pi), which is 
already known according to Invariant (14). Further, this trapezoid is the one that contains pi since 
Pi G or{pi)p. Starting at pi in this trapezoid, we move along the segment pip'i, and enter /exit 
trapezoids in HM(bay{cd)) one after another, until we encounter either p'^ or an edge of d for the 
first time. In this way, we can determine whether pip'i intersects d. Further, if pip'i intersects 
d, then the first such intersection point, denoted by z, is also found in this moving process; if 
B']^{rj,ri+i) does not intersect d, then the trapezoid of HM{bay{cd)) containing the point p'^ is 
determined. It is easy to see that the running time of the above moving procedure is proportional 
to the number of trapezoids in HM{bay{cd)) that we visit when moving along pip^. We will analyze 
the total running time of the moving process in a global manner later. 

Depending on whether i3^(rj, rj+i) (= pip'i) intersects 9, there are two cases to consider. 

If B'j^j(rj,ri^i) intersects d, then we have found the first intersection point z of B'j^j{rj,ri^i) 
and d (see Fig. 16). Note that z must be after p on d. Also, note that the Voronoi region VD(rj) 
has been computed and set as inactive, and thus pi lies on the boundary of the current active 
region of bay{cd) (see Fig. 16). Similarly as before, the line segment piz divides the current active 
region of bay{cd) into two simple polygons; one of them, say bay' , contains the point p. Then, the 
Voronoi regions of the roots that define the rays in S form a decomposition of bay', and we use 
a procedure similar to the splitting procedure discussed earlier to compute this decomposition of 
bay' , i.e., by considering the rays in S from bottom to top. However, it is possible that the target 
points of some rays in S have not been computed yet. Recall that all target points of the rays in 
S before p (= tp{pi)) have been computed. But, if the target point of a ray in S is on d{p, z), then 
it is not yet known. To compute these target points, we simply scan d{p, z) from p to z. Again, 
by Invariant (9), the vertical rays in S are target-sorted. Hence this computation can be done in 
linear time in terms of the number of edges of d{p, z) and the number of target points found during 
this process. In addition, we set the region bay' and all rays in S as inactive, and pop all rays out 
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Figure 17: Illustrating an example that the point p'l {— or{p*)) is on the bottom edge of Rec{pi,qi). 

of S {S becomes empty). Finally, we move the reference point p* to z, and continue with the next 
ray pi+i- Again, since S is empty, similar to the analysis for the subcase 2(a) when S is empty, all 
invariants of the algorithm hold. We omit the details of the proof. 

If B'j^{rj,riJ^i) {= pip'i) does not intersect d (see Fig. 17), then as shown above, B'j^j{rj,riJ^i) 
appears entirely in SPM{T) since it is contained inside bay{cd). Again, the point p'l is on either 
the right edge or the bottom edge of Rec{pi,qi) (two cases). We discuss these two cases below. 
Recall that the trapezoid of HM{bay{cd)) that contains p'^ has been computed. Recall that the ray 
at the top of S is p, lying on B{rj,ri) with i > j. 

We first discuss the case when p'^ is on the bottom edge of Rec{pi,qi) (see Fig. 17). Let p* be 
the vertical ray originating at p'^ and going south, which is on B{rj,ri^i) by Observation 2. We 
pop p out of S and push p* onto the top of S, and set p as inactive and p* as active. We move p* 
to p (= tp{pi)). We then continue to consider the next ray pi+i £ ^. Lemma 17 below shows that 
all invariants of the algorithm hold. 

Lemma 17 At the moment right before the next ray pi+i is considered, all invariants of the algo- 
rithm hold. 

Proof: Let ^ be the moment right before the next ray pi+i is considered, and ^' be the moment 
right before the ray pi is considered. Thus, the change to S from the time ^' to ^ is that p is 
popped out and p* is pushed in. Recall that at the moment all invariants of the algorithm hold. 
Our goal is to prove that all invariants still hold at the moment ^. We assume that S has at least 
two rays at the moment ^ (otherwise, all invariants related to any other rays in S than p* hold 
trivially). Let p' be the second ray from the top of S (i.e., right below p* in S) at the moment ^. 
Then p' is also the second ray from the top of S at the moment ^' (i.e., right below p in S). See 
Fig. 17 for an example. 

Invariant (1) simply follows since p* is vertically going south. 

For Invariant (2), it suffices to show that or{p*) is to the southwest of or(p'). At the moment 
since the ray p at the top of S is on B{rj,ri) with i > j, by Invariant (6), rj is to the left of 

or{p'). Since or{p*) G i?ec(rj, rj+i) is to the left of rj, we obtain that or{p*) is to the left of or{p'). 

Further, at the moment by Invariant (2), or{p) is to the southwest of or(p'). Since or{p*) is 

below or{p), or{p*^) is below or{p'). Since or{p*) is both below and to the left of or{p'), we obtain 

that or(p*) is to the southwest of or{p'). Invariant (2) thus follows. 

For Invariant (3), we need to show that or(pi^i) is to the southwest of or{p*). By Lemma 10, 

or{pi^i) is to the southwest of or{pi). Since or{p*) is to the right of or{pi), or{p*) is also to the 
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right of or(pi^i). By Lemma 8(3), or(pj_|-i) is below rj+i. Hence, or(pj_|-i) is below the rectangle 
Rec{pi, qi) and thus below or{p*). Since or(pj_|_i) is both below and to the left of or{p*), we obtain 
that or{pi+i) is to the southwest of or{p*). Thus, Invariant (3) follows. 

Invariant (4) simply follows since p* C B(rj, ri+i) and /Sj+i C B{ri+i,ri+2)- 

For Invariant (5), we need to consider both p* and Pi+i- For p*, since p* C B(ry , r,;+i), we 
need to show that the boundary portion of the Voronoi region VD{ri-^-i) (resp., VD(rj)) from Vi 
(resp., Vj) to or{p*) has been computed. For this, recall that the boundary portion of VL>(ri+i) 
between Vi and pi is a common boundary of VX)(rj_|_i) and VD(rj), which has been computed. 
We denote this boundary portion by a. Also, the boundary portion of VD{rj) between vj and 
pi has been computed; wc denote this boundary portion by /3. Further, after we find the point 
p'l, the line segment pip'i has also been obtained. Since pip[ appears entirely in SPM{T), the 
boundary portion of VD(rj+i) between Uj and or{p*) (= p[) is the concatenation of a and pip'i, 
which has been computed. Similarly, the boundary portion of VD{rj) between vj and or{p*) is the 
concatenation of /3 and Pip'i, which has been computed too. Thus, the case for p* holds. 

For the ray pi+i, which is the ray to be considered next by the algorithm, we need to show 
that the boundary portion of the Voronoi region VX'(ri_|_i) from Wj+i to or(/9i+i), which is also the 
boundary portion of the Voronoi region T/D(rj+2) from Vi+i to or{pi+i), has been computed. This 
simply follows from Observation 6. 

In summary, Invariant (5) holds. 

For Invariant (6), it suffices to show that or{p*) is to the left of rj and to the right of rj+i. Recall 
that or{p*) is on the rectangle Rec{pi,qi), qi is to the southeast of pi, and qi is the intersection 
of the vertical line through Vj and the horizontal line through rj+i. As shown above, rj+i is to the 
left of pi. Since or{p*) is to the right of pi, or{p*) is to the right of r^+i. Since Rec{pi,qi) is to the 
left of Tj, or{p1) is to the left of rj. Invariant (6) thus holds. 

For Invariant (7), we need to show that Vj is to the left of all rays in 5 \ {p*}- At the moment 
the ray p C B{rj,ri) is at the top of S with i > j; thus by Invariant (7), rj is to the left of all 
rays in 5" \ {p}. Since S \ {p} = S \ {p*}, Invariant (7) still holds at the moment ^. 

For Invariant (8), recall that p' is the second ray from the top of S at both the moments ^ and 
We assume p' lying on B{rt, r^) with t' > t. To prove Invariant (8) held at ^, it suffices to show 
t' = j since p* C -B(rj,rj+i). At the moment since p C B{rj,ri) is the ray at the top of S, by 
Invariant (8), we have j = t' . Thus, Invariant (8) still holds at the moment ^. 

For Invariant (9), it suffices to show that tp{p*) is after tp{p') on d. Intuitively this is true due 
to the following facts: There is a path inside bay{cd) from Vi to or{p*) (i.e., the concatenation of 
Vior{pi), or{pi)pi, and pip'i), and both p* and p' are vertical, and p* is to the left of p'. A detailed 
analysis is given below. 

First, it is easy to see that tp{p*) must be after the point p (= tp{pi)). The target point tp{p') 
may be after p or before p. If tp{p') is before p, then we are done. Thus, we consider the case of 
tp{p') being after p. By Invariant (2) (at the moment 4) proved above, or{p') is to the northeast 
of or{p*). Thus, the ray p' must cross pi before it hits d at tp{p'); in other words, the two line 
segments or{p')tp{p') and or{pi)tp{pi) intersect inside bay{cd). Further, since p' is to the right of 
p* and pi is to the left of p*, the intersection point of or{p')tp{p') and or{pi)tp{pi) is on pitp{pi). 

Recall that Vior{pi) is either a single point or a line segment that is in VD{bay{cd)) and does 
not intersect d. Consider the region in bay{cd) bounded by Vior{pi), or{pi)p, and d{p,d), which 
we denote hy Z. It is easy to see that Z is a simple polygon. Let a be the concatenation of pip'^ 
and p'itp{p*). Note that a is entirely inside Z except that its two endpoints are on the boundary 
of Z, i.e., pi € or{pi)p and tp{p*) G d{p,d). Thus, a divides Z into two simple polygons; one of 
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them contains pip as an edge, which is denoted by Z' . Since the intersection of or{p')tp{p') and 
or{pi)tp{pi) is on pitp{pi), the ray p' intersects Z'. By Invariant (7) (at the moment ^) proved 
above, the root rj is to the left of p'. Thus, p' cannot intersect the curve a. Hence, the target point 
tp{p') must be on the boundary of Z' fl d{p, d), which is on d{p, tp{p*)). Thus, tp{p*) is after tp{p'), 
and Invariant (9) follows. 

For Invariant (10), we need to show that if p^+i is vertical, then the target point tp{pi^i) is 
after tp{pl) on d. By Invariant (3) (at the moment ^) proved above, or(pi+i) is to the southwest 
of or{p*). Let Z be the simple polygonal region in bay{cd) bounded by Vior{pi), or{pi)pi, pip'i, 
p'ltpiPi), d{tp{p*), d), and dvi. Regardless of whether or{pi+i) = Uj+i, the origin or{pi+i) of pj+i is 
in Z since or{pi-^-i) is to the southwest of or{p*) = p'-^. Further, since both p| and pi+i are vertical, 
tp{pi-\.i) must be on d{tp{pl),d). Invariant (10) thus follows. 

For Invariant (11), it suffices to show that tp{p*) is after p* since tp{p*) has not been computed. 
Since tp{p*) is after p = tp{pi) (= p*), Invariant (11) simply follows. 

For Invariant (12), we need to show that the target point tp(/>j+i) is after p* {= p = tp{pi)). Let 
Z be the simple polygonal region in bay{cd) bounded by ViOr{pi), or{pi)p, d{p, d), and dvi. Clearly, 
or(pi+i) is in Z. Further, since or(pi+i) is to the southwest of or(pi), regardless of whether pj+i is 
vertical or horizontal, tp{pi+i) must be on d{p*,d). Thus, Invariant (12) holds. 

For Invariant (13), suppose I is the smallest index with / > i such that € ^ and pi is 
horizontal. We need to prove that tp{pi) is after p* {= p = tp{pi))- Consider the simple polygon 
Z defined above for proving Invariant (12). Since pi is horizontal, by Lemma 11, pi is below rj-i-i. 
Thus, it is easy to see that or{pi) is in Z and tp{pi) is on d{p*,d). Hence, tp{pi) is after p*, and 
Invariant (13) holds. 

For Invariant (14), we need to show that the trapezoid of HM{hay{cd)) that contains or(/9i+i) 
is known, which is true by Observation 6. 

We conclude that all invariants of the algorithm still hold at the moment ^. □ 

For the purpose of discussing the analysis of the running time of our algorithm later, we call 
the ray p| the termination vertical ray of the (horizontal) ray pi. 

We have finished the discussion for the case when p'^ is on the bottom edge of Rec{pi, qi). 

We then discuss the case when the point p[ is on the right edge of Rec(j)i,qi) (see Fig. 18). 
Denote by pu the horizontal ray originating at p'l and going east, which is on B{rj,ri+i) by 
Observation 2. Then, we pop p out of S and set p as inactive. Also, we set pn as active and move 
the reference point p* to p (= tp{pi)). Finally, we let pn be the next ray to be considered by the 
algorithm (note that pn is not in 5'). Lemma 18 below shows that all invariants of the algorithm 
hold. Recall that the trapezoid of HM{bay{cd)) that contains p'l has been computed. 

Lemma 18 At the moment right before the next ray pn is considered, all invariants of the algorithm 
hold. 

Proof: Let ^ be the moment right before the next ray pn is considered, and ^' be the moment 
right before the ray pi is considered. Thus, the only change to S from the time ^' to ^ is that p is 
popped out. At the moment all invariants of the algorithm hold. Our goal is to prove that all 
invariants still hold at the moment ^. We assume 5 7^ at the moment ^ (otherwise, all invariants 
related to any rays in S hold trivially). Let p' be the ray at the top of S at the moment ^. Then 
p' is the second ray from the top of S (i.e., right below the ray p in S) at the moment Refer to 
Fig. 18 for an example. 

Invariants (1) and (2) simply hold. 
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Figure 18: Illustrating an example that the point p'l (= or{pii)) is on the right edge of i?ec(pi, gi). 

For Invariant (3), we need to show that or{pii) is to the southwest of or{p'). By Invariant (2) 
at the moment or{p) of the ray p at the top of S is to the southwest of or{p'). Since or{pii) = p'^ 
is below or{p), or(pji) is below or{p'). Also, by Invariant (7) at the moment since the top ray p 
in S is on B{rj,ri) (with i > j), r.j is to the left of p' (which is vertical). Since p'l is on the vertical 
line through Vj, p\ = or{pii) is to the left of or{p'). Hence, or{pii) is to the southwest of or{p'), 
and Invariant (3) follows. 

For Invariant (4), suppose p' is on B{rt,rtf) with t' > t; we need to show j = t' since pn C 
B{rj,ri-^-i) is the next ray to be considered by the algorithm. At the moment ^' , p C B(rj,ri) (with 
i > j) is at the top of S and p' is the second ray from the top of S; thus, by Invariant (8) at the 
moment S,' , j = t' . Invariant (4) hence follows. 

For Invariant (5), since no new ray is pushed onto S, we only need to consider the ray pn. The 
proof is the same as that for Invariant (5) (for the ray p|) in the proof of Lemma 17, and we omit 
it. 

Invariants (6), (7), (8), (9), and (11) trivially hold since no new ray is pushed into S. 
Invariant (10) simply follows since pn is the next ray to be considered by the algorithm and pn 
is not vertical. 

For Invariant (12), we need to show that the target point tp{pii) is after p* {= p = tp{pi)). 
Consider the simple polygonal region Z in bay{cd) bounded by Vior{pi), or{pi)p, d{p,d), and dvi. 
It is easy to see that or{pii) is in Z and tp^pn) is on d{p*,d). Thus, tp^pn) is after p* . 

For Invariant (13), suppose I is the smallest index with / > i such that pi (z ^ and pi is 
horizontal. We need to prove that tp{pi) is after p* {= p = tp{pi)). Consider the simple polygon 
Z defined above for proving Invariant (12). Since pi is horizontal, by Lemma 11, pi is below ri^i. 
Thus, it is easy to see that or{pi) is in Z and tp{pi) is on d{p*,d). Hence, tp{pi) is after p* , and 
Invariant (13) holds. 

For Invariant (14), recall that the trapezoid of HM{bay{cd)) that contains p'^ (= or{pii)) has 
been computed, and thus Invariant (14) holds. 

We conclude that all invariants of the algorithm hold at the moment ^. □ 

For analysis, we refer to the ray pn as a successor horizontal ray of the (horizontal) ray pi. 

This finished the discussion for the case when p[ is on the right edge of Rec{pi,qi). 

Again, pn is the next ray to be considered by the algorithm. Although our earlier discussion 
on the algorithm processing the next ray is mostly on processing a ray pj € \I', the processing for 
Pii (^ ^) is the same, and the proof for all invariants is also very similar. In particular, there may 
also be a termination vertical ray or a successor horizontal ray generated at the end of processing 
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Figure 19: Illustrating the first two successor horizontal rays pn and pi2 of a horizontal ray £ Vl/. 

Pil, which wc stiU refer to as a termination vertical ray or a successor horizontal ray of pi. It is 
easy to see that a horizontal ray pi may lead to multiple successor horizontal rays but at most one 
termination vertical ray, i.e., a successor horizontal ray may generate another successor horizontal 
ray (e.g., see Fig. 19), but a termination vertical ray does not generate another ray. 

One might be curious about why the roles of horizontal rays and vertical rays arc quite different 
in our above algorithm, while the Li metric docs not prefer one of these two directions over the 
other. The asymmetric roles of these two directions are related to the order of pi,p2, • • • ,Pk-i in 
which we process these rays. If one uses a reversed order (i.e., pk-i,Pk-2, ■ ■ ■ , Pi) in the processing, 
then the roles of these two types of rays will be reversed. 

For the purpose of analyzing the running time of the algorithm later, we discuss more details 
related to the successor horizontal rays of a horizontal ray G We process the first successor 
horizontal ray pu of pi in the same way as pi. After pn is processed, we may obtain another 
successor horizontal ray pi2. In general, assume all successor horizontal rays we obtain for pi are 
Pil, Pi2, ■ ■ ■ , Pit, ordered by the time when they arc produced (see Fig. 19). Then, after the last ray 
Pit is processed, we may or may not obtain the termination vertical ray p^. For example, when 
processing pit, if 5 = 0, then no termination vertical ray is generated. In either case, after pit is 
processed, we continue to consider the next ray p^+i G 

Let pio = Pi. For each 1 < tt; < t, we define the points Pw+i, Qw+i, and p^+i for the ray piw 
similarly to the points pi, qi, and p'^ for piQ (sec Fig. 19). Note that when processing pit, depending 
on the specific situations, the points pt+i, qt+i, and p^_^_i may not exist (e.g., if 5 = 0). In the 
following, we assume they exist (otherwise, the analysis is actually simpler). 

It is easy to see that for each 1 < w < t, the ray Pi^w-i contains the top edge of the rectangle 
Rec{pw,p'yj) and the ray piw touches the bottom edge of Rec{pw , p'^j) . In addition, the ray pa 
contains the top edge of Rec{pt+i,p[_^_i). In other words, pio (= pi), Rec{pi,p[), pn, Rec{p2,p'2), Pi2, 
. . ., Rec{pt,p[), Pit, Rec{pt+i,p[_^_i) are ordered from high to low and left to right (see Fig. 19). Thus, 
no two different rectangles in the sequence above intersect in their interior. Actually, the rectangles 
Rec{p\,p'-y).,Rec{j)2,p'2), ■ ■ ■ ,Rec(j)tJ^i,p't^i) are ordered from northwest to southeast. Further, all 
successor horizontal rays and rectangles involved are higher than rj+i. To see this fact, note that 
for each 1 < w < t + 1, the point p'^ is higher than the point qyj and qyj is on the horizontal 
line through r^+i (see Fig. 19). Thus, all these rectangles are contained in the horizontal strip 
between the horizontal line containing pi and the horizontal line through r^+i ; we denote this strip 
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by HStrip(pi). Recall that during our algorithm, the horizontal visibility map HM{hay{cd)) is 
utilized as a guide and wc often move among its trapezoids. When computing and processing 
these successor horizontal rays, we always follow HM{bay{cd)), e.g., for each w = 0,1, . . . ,t, we 
utilize HM{hay{cd)) to compute the target point tp{piw) of the ray piw, to determine whether 
Pw+ip'w+i intersects d, and to find the trapezoid in HM{bay{cd)) that contains p'^j^i = or(/9j^^+i). 
The discussion above implies that the time for processing all successor horizontal rays of pi is 
proportional to 0{t) plus the number of trapezoids in HM{hay{cd)) that intersect the horizontal 
strip HStrip{pi) as well as the time for computing the target points of some (vertical) rays in S. 

In addition, during this process, each of the t successor horizontal rays piw of pi corresponds to 
a ray in S that is popped out. Thus, there are t vertical rays popped out of S for pi. But, at most 
one ray, i.e., the termination vertical ray pT, is pushed onto S for pi. 

We have finished the description of our algorithm for computing SPM{bay{cd)), which is sum- 
marized by the pseudo-code of Algorithm 1. 

5.3.3 The Time Complexity 

It remains to analyze the running time of the algorithm. First, we show the following lemma. 

Lemma 19 The total number of rays ever contained in the stack S throughout the entire algorithm 
is at most k. Once a ray is popped out of S, it will never be pushed back in again. 

Proof: When processing each ray pi G ^, if it is vertical, then we push it onto S; if it is horizontal, 
then as shown above, although there may be multiple successor horizontal rays of pi, at most one 
ray, i.e., the termination vertical ray, is put into S. Further, according to our algorithm, once a ray 
in S is popped out, it will never be considered again, and thus never be put into S again. □ 

We then discuss the total time for computing the target points for all vertical ray shootings in 
the entire algorithm. We use a reference point p* on d and the vertical visibility map VM{bay{cd)) 
for this purpose. To conduct the vertical ray shootings, because the rays involved are always target- 
sorted, we simply scan the edges in a portion of d between p* and another point p that is after 
p* on d. Further, when such a scanning is done, we always move p* to p. This implies that any 
portion of d is scanned at most once in the entire algorithm. In addition, the number of all vertical 
ray shootings is at most k. This is because each vertical ray involved is from S, and by Lemma 19, 
the number of rays ever contained in S is at most k. Therefore, the total time for computing the 
target points of all vertical rays in the entire algorithm is 0{n' + k). 

For each ray G if it is vertical, then processing it takes 0(1) time, i.e., pushing pi onto 
S. If it is horizontal, then assume that pi has t successor horizontal rays. We have discussed that, 
besides the procedure for computing their target points, the time for processing these t successor 
horizontal rays is proportional to t plus the number of trapezoids in HM{bay{cd)) intersecting the 
horizontal strip HStrip{pi). We have also shown that each successor horizontal ray corresponds 
to a ray in the stack S that is popped out. Since there are at most k rays ever contained in S 
by Lemma 19, the total number of successor horizontal rays in the entire algorithm is at most k. 
On the other hand, consider two different horizontal rays pi and pj in ^. We claim that the two 
horizontal strips HStip{pi) and HStrip{pj) do not intersect each other in their interior. WLOG, 
assume i < j. Indeed, the strip HStrip{pi) is above the horizontal line through the root r^+i and 
HStrip{pj) is below the ray pj. Since pj is horizontal and j > i, hy Lemma 11, pj is below rj+i. 
Our claim thus holds. The above claim implies that, besides the time for computing their target 
points, the time for processing all successor horizontal rays in the entire algorithm is proportional 
to the total number of trapezoids in HM{bay{cd)) plus k, which is 0{n' + k). 
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The algorithm performs totally 0{k) horizontal ray shootings, for computing the target points 
of the horizontal rays in ^ and their successor horizontal rays. Using HM{bay{cd)) and based on 
the fact that we already know (i.e., have computed) the trapezoid of HM{bay{cd)) containing the 
origin of each such horizontal ray, all such horizontal ray shootings can be done in 0{k) time. 

In summary, the total running time of our algorithm for computing the shortest path map for 
the bay hay{cd) is 0{n' + m!) (where m! = k — 1 \s the number of SPM{M) vertices on cd). It is 
easy to see that the size of this SPM is 0{n' + m') (e.g., since the running time is 0{n' + m')). 

Theorem 5 thus follows. 

6 Computing a Shortest Path Map for a Canal 

In this section, we show how to compute a shortest path map for a canal, which uses our shortest 
path map algorithm for a bay in Section 5 as a main procedure. 

Consider a canal canal{x, y) with x and y as the corridor path terminals and two gates xd and 
IJz (e.g., see Fig. 6). There may be SPM{Ai) vertices on both gates. Let mi (resp., ni2) be the 
number of SPM{M) vertices on xd (resp., yz), and n' be the number of obstacle vertices of the 
canal. We show that a shortest path map for the canal can be computed in 0{mi + m2 + n') time. 
Let R\ (resp., R2) be the set of roots whose cells in SPM{M.) intersect xd (resp., yz). 

Recall that we have defined wavefront incoming/ outgoing terminals in Section 4.2. Namely, 
consider the corridor path terminals x and y of canal{x,y). It is possible that y has a shortest 
path from s via x (i.e., this path contains the corridor path of canal(x,y)), in which case there is 
a "pseudo-cell" in SPM{M.) with x as the root and y being the only other point in this "pseudo- 
cell" ; then x is a wavefront incoming terminal and y is the wavefront-outgoing terminal. If neither 
y has a shortest path from s via x nor x has a shortest path from s via y, then both x and y are 
wavefront-incoming terminals. In this case, there is a point on the corridor path of canal{x, y) that 
has two shortest paths from s, one via x and the other via y (we will use this property to compute 
an SPM for canal{x,y)). 

Note that for the two terminals x and y, either both of them are wavefront-incoming terminals, 
or only one of them is an wavefront-incoming terminal and the other is an wavefront-outgoing 
terminal. Below, we first discuss the former case; the algorithm for the latter case is very similar. 

6.1 Both X and y are Wavefront-incoming Terminals 

If both X and y are wavefront-incoming terminals, by the properties of the corridor path, there is a 
point p* on the corridor path of canal{x,y) such that there exist two shortest paths 7ri(s,p*) and 
712(3, p*) from s to p* with x G 7ri(s,p*) and y G 7r2(s,p*). The point p* can be found in 0{n!) time 
since we know the shortest path distances from s to x and to y. 

Let VD (canal [x, y), Ri) be the (additively) weighted Voronoi diagram of canal(x,y) with re- 
spect to the root set Ri, i.e., we treat canal(x, y) as a bay with the gate xd. As defined in Section 
5, VD(canal(x, y), Ri) is the Voronoi decomposition of canal(x, y) with respect to the roots in Ri. 
Similarly, let VD (canal (x, y), R2) be the weighted Voronoi diagram of canal(x,y) with respect to 
the root set i?2- Using our algorithm in Section 5, VD (canal (x, y), Ri) and VD(canal(x,y), R2) 
can be computed in totally 0(mi + m2 + n') time. Denote by VD(canal(x, y),Ri,R2) the weighted 
Voronoi diagram of canal (x,y) with respect to the roots in Ri L) R2. As shown in Section 5, af- 
ter VD (canal (x, y), Ri, R2) is computed, an SPM on canal(x,y) with the source s can be built 
in 0(mi + 777,2 + n') time. Thus, the key is to compute VD (canal (x, y), Ri, R2). Below, we show 
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how to compute VD (canal (x, y), Ri, R2) in 0(mi + m2 + n') time with the help of the point p*, 
VD {canal {x, y), Ri), and VD (canal (x, y), R2)- 

To compute VD(canal(x,y), Ri, R2), our strategy is to find a "dividing curve" in canal(x,y) 
that divides canal{x,y) into two simple polygons Ci and C2, such that each point in Ci has a 
shortest path from s via a root in Ri and each point in C2 has a shortest path from s via a root in 
i?2- Further, each point on the dividing curve has two shortest paths from s, one path containing 
a root in Ri and the other path containing a root in R2. After finding Ci and C2, we simply apply 
the algorithm in Section 5 on Ci and Ri to compute the weighted Voronoi diagram of Ci with 
respect to Ri, i.e., VD(Ci, Ri). We similarly compute VD(C2,i?2)- Then, VD (canal (x, y), Ri, R2) 
consists of VD(Ci, Ri) and VL'(C2, -^2)- Thus, our remaining task is to compute a dividing curve 
in canal(x,y), which wc denote by 7. 

Note that the point p* G 7. Computing 7 can be done in 0(n' + mi + m2) time by a procedure 
similar to the merge procedure of the divide-and-conquer algorithm for computing the Voronoi 
diagram of a set of points in the plane [34]. The details are given below. 

To compute 7, we start at the point p* and trace 7 out by traversing some corresponding cells 
in VD(canal(x,y), Ri) and in VD (canal (x, y), R2) simultaneously. Specifically, we first compute a 
triangulation of VD (canal (x, y), Ri), denoted by Trii, and a triangulation of VD (canal (x, y), R2), 
denoted by Tri2 (this can be done in linear time [3] since each cell of VD (canal (x, y), Ri) and 
VD (canal (x, y), R2) is a simple polygon). Since p* is in a triangle (say, trii) of Trii and is in a 
triangle (say, tri2) of Tri2, we find trii in Trii and tri2 inTri2- From the cell of VX'(cana/(x, y), 
(resp., VD(canal(x, y), R2)) that contains trii (resp., tri2)-, we obtain the root ri (resp., r2) of that 
cell. We then move along the bisector B(r\,r2) inside canal(x,y), starting at p* and going in each 
of the two directions along B(ri,r2)- As following a line segment or a ray of B(ri,r2) in a direction, 
we determine, in 0(1) time, which of trii or tri2 that we exit first. As we cross from one triangle 
tri (say, in Trii) to the next triangle tri', we check which of the following cases occurs: (i) The 
next triangle tri' (in Trii) is contained in the same cell of VD (canal (x, y), Ri) as that containing 
tri; (ii) tri' is contained in a different cell of VD (canal (x, y), Ri) than that containing tri; (iii) the 
movement touches the boundary of canal(x,y) (thus tri' does not exist). In Case (i), we continue 
to follow the same bisector (say, B(ri,r2))- In Case (ii), we find the root (say, r'^) of the next cell of 
VD(canal(x, y), Ri); then we compute a new bisector (say, B(r'i,r2)), and our movement continues 
along B(r'i,r2)- In Case (iii), the movement reaches an end of 7 (on the boundary of canal(x,y)). 
The dividing curve 7 is the concatenation of the portions of the bisectors thus traversed. 

Due to the properties of the cells of VD(canal(x, y),Ri) and VD(canal(x, y), R2), our movement 
above can visit each triangle of Trii and Tri2 at most once, taking 0(1) time per triangle visited. 
Thus, the partition curve 7 is computed in 0(n' + mi + 777.2) time. 

In summary, in this case, an SPM on canal(x, y) can be computed in 0(n' + mi + 7712) time. 

6.2 Only One of x and 7/ is a Wavefront-Incoming Terminal 

In this case, exactly one of x and y is a wavefront-incoming terminal. The algorithm is similar to 
that for the former case. The only difference is on how to find a point p* on the dividing curve 7 
because in this case no such a point p* can be on the corridor path of canal(x, y). 

WLOG, we assume that x is a wavefront-incoming terminal and y is not. Then each point on 
the corridor path (including y) has a shortest path from s via x. Further, the shortest path through 
X passes y and goes to the outside of canal(x,y), which means that y is the root of a cell C(y) in 
SPM(M). If the canal gate yz is completely contained in the cell C(y), then it is easy to see that 
VD(canal(x, y),Ri) is VD(canal(x, y), i?i, i?2)- Otherwise, as in the former case, we need to find a 
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dividing curve 7 to divide canal{x, y) into two polygons Ci and C2 such that each point in Ci has 
a shortest path from s via a root in R\ and each point in C2 has a shortest path from s via a root 
in R2. To obtain 7, the key is to find a point p* G 7. Since the canal gate 'yz is not completely 
contained in C(y), there must be a point q on yz that is on the common boundary of C{y) and 
another cell C(r) in SPM{M). We claim that q is on 7. Indeed, note that r is in i?2- Hence there 
is a shortest path TTi{s,q) from s to g that contains x, the corridor path in canal{x,y), and the 
line segment yq, and there is another shortest path 7r2{s,q) from s to g via the root r G i?2- In 
other words, q has two shortest paths from s, one via a root in Ri and the other via a root in R2. 
Therefore, q is on the dividing curve 7. The rest of the algorithm is similar to that for the former 
case. 

In summary, in this case, an SPM on canal{x, y) can also be built in 0{n' + mi + 771-2) time. 
Therefore, a shortest path map SPM on canal{x, y) can be computed in 0{n' + mi + 777-2) time. 
Similarly, the size of this SPM is 0{n' + m\ + m2)- 
Theorem 6 thus follows. 

7 Applications of Our Shortest Path Algorithms 

In this section, we extend our techniques to solve some other problems. 
7.1 The Li Geodesic Voronoi Diagram 

Given a set V oi h polygonal obstacles of totally n vertices and a set of m point sites, the Li-GVD 
problem aims to construct the Li geodesic Voronoi diagram of for the m point sites. Denote by 
GVD{V) the Voronoi diagram that we want to construct. 

Mitchell's algorithm [29, 30] can be modified to compute GVD{V) in 0{{n + 7n)log(77 + m)) 
time. Namely, instead of initiating a wavelet at a single source, the modified algorithm for GVD('P) 
initiates a wavelet at each point site. The rest of the algorithm remains the same as before. 

We can also extend our SPM algorithm in a similar way to compute GVD{V). Generally, since 
our algorithm makes use of Mitchell's algorithm [29, 30] as a main procedure when computing the 
shortest path map SPM{M.) for the ocean A^, to compute C!VD{V), we can simply replace Mitchell's 
algorithm by its modified version for computing Li geodesic Voronoi diagrams. More specifically, 
our algorithm for computing GVD{V) has the following steps. (1) Compute a triangulation of the 
free space, in which the 777- point sites are treated as 777- point obstacles. (2) Compute the corridor 
structure on V and the m point obstacles that consists of 0{m + h) corridors, which partition the 
plane into a set V' of 0{m + h) convex polygons of totally 0(77 + 777-) vertices. (3) Compute the 
core set core{V') for the convex polygons in V' . (4) Apply Mitchell's modified algorithm [29, 30] 
to compute the Li geodesic Voronoi diagram GVD{core(V')) on the core set core{V'). (5) Based 
on GVD{core{V')), compute the Li geodesic Voronoi diagram GVD{V') on the convex polygon set 
V' . Although we have multiple sources, this step is the same as before (i.e., as in Lemma 5). (6) 
Based on GVD{'P'), compute the Voronoi regions in all bays and canals, as in Sections 5 and 6. 
Again, the algorithms for this step are as before, i.e., as the algorithms in Sections 5 and 6. We 
then obtain the final Li geodesic Voronoi diagram GVDiV). 

To analyze the running time. Steps (1), (2), and (3) are the same as before except that the 
number of obstacles becomes m + h. Specifically, the triangulation in Step (1) takes 0(77 + {h + 
777) \og^'^''{h + 77i)) time [2]. Steps (2) and (3) together take 0(77 + {h + m) log(/i + 777)) time. Step 
(4) takes 0{{m + h)log{m + h)) time since the core set coreiT") has totally 0(777 + h) vertices. 
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Steps (5) and (6) arc also the same as before, which take hnear time, i.e., 0{n + m). Therefore, 
the entire algorithm takes 0{n + [h + m) \og^^''{h + m)) time, which is dominated by the time of 
the triangulation procedure in Step (1). 

Theorem 8 The Li geodesic Voronoi diagram of m point sites among a set of h pairwise disjoint 
polygonal obstacles of totally n vertices in the plane can be computed in 0{n+{h+m) log^^'^(/i+m)) 
time (or 0{n + {h + m) \og{h + m)) time if a triangulation is given). 

If the m point sites are all inside a simple polygon, then Theorem 8 leads to the following result. 

Corollary 2 The Li geodesic Voronoi diagram of a set of m point sites in a simple polygon can 
be computed in 0{n + mlog^'^^'^m) time (or 0{n + mlogm) time if a triangulation is given). 

Note that the currently fastest known GVD algorithm for the Euclidean version of the single 
simple polygon case runs in 0((n + m) log(n + m)) time [32]. 

Remark. Since the given m sites are points, there is an alternative triangulation algorithm 
that may be faster (than simply applying the algorithm in [2]) in some situations. The algorithm 
works as follows: (1) Compute the triangulation of the free space without considering the m sites; 
(2) find the triangles in the triangulation that contain those m sites (e.g., by a point location data 
structure); (3) triangulate those triangles that contain at least one point site by considering the 
point sites as obstacles. It is easy to see that this algorithm takes 0(n + m log n) time in the single 
polygon case and 0{n+h log"^^^ h+m log n) time in the polygonal domain case. Therefore, using this 
triangulation algorithm, the geodesic Voronoi diagram can be constructed in 0(n+m(logn+logm)) 
time in the single polygon case and in 0{n + h log^"*"^ h + m log n+ {h + m) \og{h + m)) time in the 
polygonal domain case. 

7.2 Shortest Paths with Fixed Orientations and Approximate Euclidean Short- 
est Paths 

As in [29, 30], our algorithms can be generalized to solving the C-oriented shortest path problem 
[37]. A C-oriented path is a polygonal path with each edge parallel to one of a given set C of 
fixed orientations. A shortest C-oriented path between two points is a C-oriented path with the 
minimum Euclidean distance. Rectilinear paths are a special case of this problem with two fixed 
orientations of and 7r/2. Let c = |C|. Mitchell's algorithm [29, 30] can compute a shortest C- 
oriented path in O (era log n) time and 0{cn) space among h pairwise disjoint polygons of totally n 
vertices in the plane. Similarly, our algorithms also work for this problem, as follows. 

We first consider the convex case (i.e., all polygons are convex). We compute a core for each 
convex polygon based on the orientations in C. Note that in this case, a core has 0(c) vertices. 
Thus, we obtain a core set of totally 0{ch) vertices. We then apply Mitchell's algorithm for the 
fixed orientations of C on the core set to compute a shortest path avoiding the cores in 0{(?h\og ch) 
time and 0{c^h) space, after which we find a shortest path avoiding the input polygons in additional 
0{n) time as in Lemma 4. Thus, a shortest path can be found in totally 0{n + (?h\ogch) time 
and 0{n + (P'h) space. For the general case when the polygons need not be convex, the algorithm 
scheme is similar to our L\ algorithm in Section 4. In summary, we have the following result. 

Theorem 9 Given a set C of orientations and a set of h pairwise disjoint polygonal obstacles of 
totally n vertices in the plane, we can compute a C-oriented shortest s-t path in the free space 
in 0{n + /ilog^^*^ h + c^hlogch) time (or 0{n + c^hlogch) time if a triangulation is given) and 
0{n-\- (?h) space, where c= \C\. 
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This also yields an approximation algorithm for computing a Euclidean shortest path between 
two points among polygonal obstacles. Since the Euclidean metric can be approximated within an 

accuracy of O if we use c equally spaced orientations, as in [29, 30], Theorem 9 leads to an 
algorithm that computes a path guaranteed to have a length within a factor + of the Euclidean 
shortest path length, where c is chosen such that 5 = 0(l/c^). 

Corollary 3 A 5 -optimal Euclidean shortest path between two points among h pairwise disjoint 
polygons of totally n vertices in the plane can be computed in 0{n + hlog^^'^ h + {l/d)hlog^) time 

(or 0{n + (l/(5)/ilog time if a triangulation is given) and 0{n + {l/5)h) space. 

8 Conclusions 

We present new algorithms for solving Li shortest path problems in polygonal domains. Our 
algorithms are optimal if the triangulation for the free space can be done optimally (i.e., T = 
0{n + /ilog/i)). In fact, our results show that building an Li shortest path map is equivalent to 
the triangulation in terms of the running time. 

Some of our techniques may be helpful on solving the Euclidean version of the problem. For 
the Euclidean version, as the Li version, a long-standing open problem is to compute a shortest 
path in 0{n + hlogh) time and 0(n) space. Hershberger and Suri [19] built an SPM of size 0(n) 
in 0(n log n) time and O(nlogn) space. Recently, Inkulu et al. announced an algorithm that 
can find an Euclidean shortest path in 0(n + hloghlogn) time [22]; we give an algorithm for the 
problem that runs in 0{n + hlog^'^'' h + k) time [8], where k is a parameter sensitive to the input 
and is bounded by 0{h'^). Note that our algorithm [8] particularly works for obstacles that have 
curved boundaries. To generalize the techniques given in this paper to the Euclidean version, some 
difficulty appears. For example, the idea of using cores does not seem to work (e.g., Lemma 4 is 
not applicable to the Euclidean version). A possible direction for solving the Euclidean version is 
to first solve the convex case with the performance desired by the open problem. If this is possible, 
then by generalizing the techniques given in this paper, it is very likely that the general case may 
also be solved accordingly, and thus the open problem can be settled although we may still have to 
suffer the 0(n + ^log^"^^ h) triangulation time. 
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